esi-rc. 


ESD-TR-69-J89 


esd  Acde'yyiuN  us- 

- — 


ESTI  Call  No 
Copy  No.  _ 


cys. 


THREE-DIMENSIONAL  CURVES  AND  SURFACES 
FOR  RAPID  COMPUTER  DISPLAY 


Ad*  !»%  /7C. 


Theodore  M.  P.  Lee 


ESD  RECORD  COPY 


30  April  1969 


SCIENTIFIC  I  mmm  MffOfrMATION  DIVISION 

(€ST-H,  WHAING  1211 


DIRECTORATE  OF  PLANNING  AND  TECHNOLOGY 

ELECTRONIC  SYSTEMS  DIVISION 

AIR  FORCE  SYSTEMS  COMMAND 

UNITED  STATES  AIR  FORCE 

L.  G.  Hanscom  Field,  Bedford,  Massachusetts 


Sponsored  by:  Advanced  Research  Projects  Agency 
Washington,  D.  C. 

ARPA  Order  No.  952 


This  document  has  been 
approved  for  public  release  and 
sale;  its  distribution  is 
unlimited. _ 


(Prepared  under  Contract  No.  FI9628-68-C-0379  by  Harvard  University, 
Cambridge,  Massachusetts.) 


LEGAL  NOTICE 


When  U.S.  Government  drawings,  specifications  or  other  data  are  used  for  any 
purpose  other  than  a  definitely  related  government  procurement  operation,  the 
government  thereby  incurs  no  responsibility  nor  any  obligation  whatsoever;  and 
the  fact  that  the  government  may  have  formulated,  furnished,  or  in  any  way  sup¬ 
plied  the  said  drawings,  specifications,  or  other  data  is  not  to  be  regarded  by 
implication  or  otherwise  as  in  any  manner  licensing  the  holder  or  any  other  person 
or  conveying  any  rights  or  permission  to  manufacture,  use,  or  sell  any  patented 
invention  that  may  in  any  way  be  related  thereto. 


OTHER  NOTICES 


Do  not  return  this  copy.  Retain  or  destroy. 


THREE-DIMENSIONAL  CURVES  AND  SURFACES 
FOR  RAPID  COMPUTER  DISPLAY 


Theodore  M.  P.  Lee 


30  April  1969 


DIRECTORATE  OF  PLANNING  AND  TECHNOLOGY 

ELECTRONIC  SYSTEMS  DIVISION 

AIR  FORCE  SYSTEMS  COMMAND 

UNITED  STATES  AIR  FORCE 

L.  G.  Hanscom  Field,  Bedford,  Massachusetts 


Sponsored  by:  Advanced  Research  Projects  Agency 
Washington,  D.  C. 

ARPA  Order  No.  952 


This  document  has  been 
approved  for  public  release  and 
sale;  its  distribution  Is 
unlimited. 


(Prepared  under  Contract  No.  FI9628-68-C-0379  by  Harvard  University, 
Cambridge,  Massachusetts.) 


FOREWORD 


This  report  describes  work  accomplished  under  Contract  F-19628- 
68-C-0379  from  June  1968  through  April  1969.  This  contract  is 
concerned  with  research  on  computer  graphics  and  computer  network¬ 
ing.  In  particular  it  is  directed  to  development  of  new  insights  into, 
the  creation,  analysis  and  presentation  of  information.  This  report 
presents  an  approach  allowing  three-dimensional  information  to  become 
an  element  in  the  repertoire  of  the  computing  machine  by  analyzing  a 
certain  class  of  three-dimensional  surfaces. 

The  report  is  based  on  a  thesis  submitted  on  April  30,  1969  by 
Theodore  M.  P.  Lee  to  Harvard  University,  Division  of  Engineering 
and  Applied  Physics,  in  partial  fulfillment  of  the  requirements  for  the 
degree  of  Doctor  of  Philosophy. 

Professor  Anthony  G.  Oettinger  was  the  principal  investigator  for 
the  contract.  Dr.  Lawrence  Roberts  was  the  ARPA  director.  Dr.  Sylvia 
R.  Mayer  was  the  ARPA  agent  at  Electronic  Systems  Division.  Lt  John 
P.  McLean  of  Electronic  Systems  Division  has  provided  technical 
guidance. 

This  technical  report  has  been  reviewed  and  is  approved. 


WILLIAM  F.  HEISLER,  Colonel,  USAF 
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ABSTRACT 


Rational  parametric  polynomial  functions  of  second  degree  or 
higher  provide  a  class  of  curves  including  all  conic  sections.  They  can 
be  generated  by  an  iterative  process  easily  implemented  in  software  or 
hardware.  The  numerical  accuracy  of  the  process  is  analyzed. 
Algorithms  for  the  specification,  display,  and  modification  of  the  curve 
are  presented.  Such  curves  are  represented  in  a  homogeneous  coordi¬ 
nate  formulation  convenient  for  computer  applications.  Three- 
dimensional  surfaces  composed  of  such  curves  are  similarly  convenient 
to  use.  Without  recourse  to  trigonometric  functions  such  classical 
surfaces  as  spheres  and  toroids  can  be  readily  described.  The  ease 
with  which  translation,  rotation  and  projective  transformations  can  be 
applied  is  exhibited.  In  particular,  we  do  not  perform  such  transfor¬ 
mations  on  the  points  of  the  surface  to  be  displayed  --  upwards  of 
several  thousand  --  but  rather  upon  the  rather  small  set  of  numbers  in 
a  4  X  4  X  4  tensor  that  represents  the  surface.  These  surfaces  are 
intended  to  be  used  in  an  interactive,  freeform  computer-aided  design 
system.  In  this  direction  we  discuss  the  enforcing  of  continuity  con¬ 
ditions  and  possible  data  structures  for  representing  the  surfaces. 
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SECTION  I 


INTRODUCTION 


1.1  GENERAL  SUMMARY 

This  dissertation  is  concerned  with  the  mathematics  and  some  of 
the  pragmatics  of  the  specification,  generation,  and  display  of  a  certain 
class  of  three-dimensional  surfaces  and  of  the  curves  out  of  which  the 
surfaces  are  formed.  It  is  assumed  that  the  reader  has  a  good  back¬ 
ground  in  linear  algebra  and  a  fair  background  in  the  art  of  computer 
graphics.  After  providing  a  survey  of  the  necessary  material  in  the 
field  of  projective  geometry  we  analyze  the  class  of  curves  of  interest 
followed  by  a  detailed  inspection  of  the  surfaces.  The  final  chapters 
present  experimental  results  as  well  as  a  proposal  for  an  evolving  data 
structure  relevant  to  the  matter  of  surface  specification. 

Recent  advances  in  computer  display  hardware  have  been  concerned 
with  achieving  the  ability  to  rapidly  draw  curves  in  two  or  three 
dimensions.  Rational  parametric  polynomial  functions  of  second  degree 
or  higher  provide  a  class  of  curves  including  all  conic  sections.  Such 
curves  can  be  generated  by  an  iterative  process  easily  implemented  in 
hardware  or  software.  The  mathematics  for  specifying,  transforming 
and  manipulating  the  curve  are  formulated  in  matrix  algebra.  In 
particular,  transformations  exist  which  alter  the  rate  of  display  of  the 
curve  without  changing  its  shape.  Cur  ves  may  be  sectioned  into  smaller 
portions  or  expanded  into  larger  continuations.  The  specific  properties 
of  rational  parametric  cubic  polynomials  are  discussed  in  detail, 
including  a  statistical  analysis  of  the  propagation  of  truncation  errors 
in  their  computation. 

The  thesis  describes  a  class  of  three-dimensional  surfaces  -- 
formed  by  bi-cubic  rational  polynomials  --  well-suited  to  computer 
display.  The  mathematical  description  of  the  surface  is  a  homogeneous 
coordinate  formulation  of  Coons'  sur  faces  providing  sufficient  degrees 
of  freedom  to  make  them  attractive  for  on-line  shape  design.  Without 
recourse  to  trigonometric  functions  such  classical  surfaces  as  spheres 
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and  toroids  are  readily  described. 

The  use  of  iterative  computation  techniques  for  generating  the 
surface  to  provide  rapid  response  is  discussed.  The  properties  of  the 
boundary  curves  and  interior  contour  curves  (all  rational  parametric 
cubics)  are  examined  The  ease  with  which  translation,  rotation  and 
projective  transformations  can  be  applied  is  exhibited.  In  particular, 
we  do  not  perform  any  such  transformations  on  the  points  of  the  surface 
to  be  displayed  --  upwards  of  several  thousand  --  but  rather  upon  the 
rather  small  set  of  numbers  in  a  tensor  that  represents  the  surface. 

Simple  continuity  constraints  and  other  ways  of  specifying  the  surface 
from  geometrical  considerations  are  studied. 

The  above  material  --  as  presented  in  Sections  II  and  III  --  is 
almost  entirely  original,  although  of  course  parallel  work  has  been 
going  on  at  a  number  of  places.  In  particular,  close  contact  has  been 
maintained  with  S.  A.  Coons,  so  it  is  difficult  to  trace  the  origin  of 
several  ideas  or  the  particular  form  of  a  derivation.  Section  I  is 
entirely  review  and  as  such  does  not  contribute  to  the  state  of  the  art. 

In  hoping  to  present  this  material  for  use  in  a  viable  interactive 
system.  Appendix  II  proposes  a  data  structure  in  which  various 
intuitive  as  well  as  mathematical  properties  of  the  surfaces  can  be 
modelled.  In  applying  some  of  the  techniques  of  artificial  intelligence 
projects  as  expressed  in  deductive  systems  it  hopes  to  free  the  user 
from  some  of  the  restrictions  of  an  algorithmic,  pre-programmed 
environment.  I  claim  originality  for  the  particular  use  of  these  con¬ 
cepts  and  for  the  application  of  an  associative  data  structure  to  the 
properties  of  surfaces,  but  must  express  the  warning  that  much  of  this 
material  is  based  on  a  cursory  exploration  of  the  possibilities  and  has 
no  claim  to  finality. 

The  basic  algorithms  for  the  display  of  curves  and  surfaces  have 
been  programmed  and  are  illustrated  throughout  the  thesis.  Most  of  the 
advanced  theory  --  such  as  specifying  a  surface  from  interior  boundary 
curves  --  has  not  been  tested  for  lack  of  a  suitably  powerful  and  con¬ 
venient  interactive  graphics  facility  with  the  capacity  of  three-dimensional 
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input.  The  two  to  three  second  response  time  for  the  computation  of  a 
Surface  on  the  available  equipment  is  just  not  adequate  for  any  future 
experiments  which  would  depend  on  the  ability  to  obtain  immediate 
feedback  for  the  manipulation  of  a  curve  or  surface  in  three  dimensions. 

Computer  graphics  --  the  technology  of  composing  and  communi¬ 
cating  two-dimensional,  primarily  non-textual,  information  with  the  aid 
of  a  computing  machine  [22]  --  is  new,  but  not  new  enough  to  let  that  be 
an  excuse  for  inefficient  techniques  or  imprecise  formulation  of  ideas. 
This  thesis  proposes  a  mathematical  formulation  for  a  certain  class  of 
curves  and  surfaces  in  three  dimensions.  This  formulation  is  suitable 
for  efficient  two-dimensional  perspective  display.  The  class  of  objects 
covered  is  sufficiently  large  as  to  be  useful  in  the  solution  of  a  non¬ 
trivial  set  of  the  problems  faced  by  people  attempting  to  use  computer 
graphics  as  an  aid  in  the  design  of  space  forms.  The  class  of  surfaces, 
in  yielding  to  direct  mathematical  characterization  and  analysis,  is  also 
interesting  in  itself.  It  is  intended  that  this  formulation  be  used  in  an 
interactive  environment  concerned  with  achieving  rapid  response  in  the 
realistic  display  of  three-dimensional  surfaces,  without  the  large  expendi 
ture  of  computing  resources  needed  in  previously  available  approaches. 
The  surfaces  are  either  to  be  used  as  a  rich  vocabulary  for  some  other 
problem  or  to  be  created  directly  as  the  objects  of  interest.  Although 
the  main  goal  of  this  work  is  the  mathematical  description  of  the 
surfaces  and  the  presentation  of  the  basic  algorithms  for  their  display, 
sufficient  practical  characterizations  are  presented  so  that  this  thesis 
may  serve  in  some  sense  as  a  reference  collection  of  directly  applicable 
techniques.  The  mathematics  is  presented  in  a  tensor  notation  easily 
transcribed  to  a  working  computer  program.  To  aid  the  problem  of 
interactively  specifying  a  surface,  characterizations  have  been  developed 
which  permit  a  small  set  of  intuitively  meaningful  geometrical  conditions 
such  as  tangent  vector  continuity  --  to  be  used  to  describe  the  desired 
surface. 

In  summary,  we  have  taken  a  simple  concise  idea,  explored  in 
detail  its  practical  implementations,  using  any  related  theory  as  an  aid  to 
understanding,  particularly  through  a  set  of  alternative  characterizations 
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useful  in  different  specific  contexts. 

Several  existing  techniques  for  modelling  three-dimensional 
surfaces  have  used  a  coordinate  dependent  representation,  such  as 
z  =  f(x,y)  where  f  is  of  some  predetermined  type,  say  formed  from 
polynomials  in  the  two  variables  x  and  y.  This  approach  may  be 
satisfactory  in  the  modelling  of  objects  with  well-defined  dependent  and 
independent  axes,  but  it  is  completely  contrived  when  applied  to  more 
general  space  forms,  such  as  an  automobile  fender  or  bumper.  For 
this  reason,  and  for  other  reasons  of  historical  accident,  we  have 
chosen  a  representation  in  which  each  coordinate  is  an  independent 
function  of  two  arguments,  the  two  degrees  of  freedom  on  a  surface. 
The  functions  have  been  chosen  to  be  computationally  efficient  -- 
polynomials  of  degree  three  or  less.  Notice  that  the  equations  of  some 
space  curves  lying  on  a  surface  can  be  obtained  by  treating  one  of  the 
arguments  as  a  constant  and  by  letting  the  other  be  the  parameter  of 
the  one  degree  of  freedom  on  the  curve. 


Instead  of  the  standard  three-dimensional  Cartesian  represen¬ 
tation  of  Euclidean  space  we  have  chosen  to  use  a  four-dimensional 
projective  representation.  This  representation,  in  its  abstract  form, 
is  taken  from  the  realm  of  projective  geometry  in  which  it  is  the  ratios 
among  components  rather  than  the  components  themselves  that 
determine  a  point  in  space.  In  this  four-dimensional  space  we  repre¬ 
sent  a  surface  (strictly  speaking,  a  part  of  a  surface,  called  a  patch) 
by  four  independent  bi-cubic  polynomial  functions,  one  for  each  of  the 
homogeneous  coordinates.  The  parameters  of  the  function  range  over 
the  unit  square,  0  |  u  g  1,  0  ^  v  ^  1.  The  four  functions  are  represented 
by  a  64  element  tensor  T  where  a  component  T  .,  is  the  coefficient  of 
the  term  u  vJ  in  the  polynomial  for  the  k’th  coordinate. 


Such  a  homogeneous  coordinate  representation  --  so  named 
because  abstractly  four  coordinates  are  used  to  represent  three  coordi¬ 
nates  with  no  preferred  mapping  from  the  one  space  to  the  other  --  has 
three  particular  advantages  over  the  Cartesian  representation  for  our 
work  in  computer  graphics.  Of  primary  importance,  the  use  of  four 
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coordinates  adds  very  convenient  extra  degrees  of  freedom  allowing 
portions  of  such  surfaces  as  spheres  or  toroids  or  of  curves  including 
all  conic  sections  to  be  represented  exactly  with  no  more  complicated 
functions  than  polynomials  of  second  degree.  In  particular,  the  extra 
freedom  provided  by  introducing  homogeneous  coordinates  allows 
surfaces  (and  the  curves  of  which  they  are  composed)  to  be  constrained 
with  respect  to  tangent  direction  and  magnitude  at  their  boundaries, 
without  sacrificing  any  internal  structure.  This  freedom  has  not  been 
possible  in  previous  related  work,  such  as  cubic  splines  or  Coons' 
bi-cubic  surfaces  [7]. 

Secondly,  placing  the  abstract  representation  of  a  surface  in 
homogeneous  coordinates  makes  the  application  of  the  proper  perspective 
projection  strictly  analogous  to  the  application  of  the  transformations  of 
rotation  and  translation.  With  the  particular  form  of  a  surface  function 
chosen,  we  thus  have  a  representation  for  objects  independent  of  the 
manner  in  which  they  are  to  be  viewed  or  composed  into  more  complex 
objects.  For  instance,  we  display  a  visual  representation  of  a  surface 
by  a  set  of  parametrically  orthogonal  curves  lying  on  the  surface.  The 
particular  set  used  and  the  coarseness  of  the  computer  representation 
of  the  curves  can  be  chosen  independently  of  the  abstract  representation 
of  the  surface,  as  can  the  particular  orientation  and  perspective  view¬ 
point,  all  such  choices  being  performed  through  simple  matrix  trans¬ 
formations  of  the  abstract  representation. 

Lastly,  the  homogeneous  nature  of  all  the  quantities  (vectors, 
matrices,  or  tensors)  used  in  a  computer  implementation  of  this  work 
allows  fixed  point  arithmetic  to  be  used  to  a  great  speed  advantage, 
with  questions  of  scaling  being  relegated  to  the  implicit  scale  factor 
present  in  all  homogeneous  arrays.  For  example,  if  the  product  of  two 
matrices  leads  to  arithmetic  overflow,  computation  of  the  product  is 
repeated  with  each  element  of  one  of  the  matrices  being  divided  by 
successively  larger  powers  of  two  until  no  overflow  occurs;  initially 
the  matrices  are  separately  normalized  such  that  the  largest  magnitude 
of  any  of  the  numbers  in  each  matrix  is  one. 
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This  thesis  is  organized  into  three  major  sections  --  Section  II 
on  curves.  Section  III  on  surfaces  and  Section  IV  discussing  ways  in 
which  the  results  of  the  previous  sections  have  been  used  in  a  practi¬ 
cal  implementation.  The  material  on  curves  does  not  depend  on  the 
material  on  surfaces,  while  the  latter  uses  explicitly  only  a  few  of  the 
major  results  of  the  former.  In  both  these  sections  the  attempt  has 
been  to  accompany  the  full  formal  presentation  with  a  few  special  case 
examples  to  indicate  the  direction  the  exposition  is  taking.  Throughout 
all  the  material  there  is  a  free  alternation  between  theoretical  con¬ 
sideration  and  explicit  practical  derivations.  Wherever  it  seems 
appropriate,  the  mathematics  has  been  presented  in  a  more  general 
form  than  would  be  strictly  necessary  for  the  practical  derivations  it 
accompanies. 

Appendix  II  proposes  a  very  general  data  structure  particularly 
appropriate  to  the  problems  of  representing  these  curves,  surfaces, 
and  the  relations  between  them.  Containing  the  least  specific  material 
in  the  thesis,  it  however  is  of  importance  in  attempting  to  represent 
the  types  of  interdependencies  in  the  various  topics,  doing  so  without 
the  necessary  but  confusing  mathematical  detail  presented  in  Sections  II 
and  III.  In  some  sense,  it  gives  an  answer  to  the  question  of  what  is  the 
purpose  of  this  research  by  presenting  a  framework  in  which  it  is 
meaningful  for  a  designer  to  ask  questions  about  the  surfaces  themselves, 
the  answers  to  which  can  be  found  from  the  mathematical  representations 
derived  so  tediously  in  the  earlier  sections. 

To  make  this  work  moderately  self-contained,  there  are  a  few 
subsections  which  are  not  strictly  speaking  original,  but  seem  of  such 
direct  relevance  as  to  merit  inclusion.  The  subsection  on  homogeneous 
coordinates  and  transformations  contains  no  new  material,  but  serves 
to  establish  a  common  background  for  the  rest  of  the  thesis.  The  basic 
ideas  for  the  iterative  generation  of  curves  are  the  joint  efforts  of 
S.  A.  Coons,  D.  Cohen,  and  the  author;  the  backwards  difference  scheme 
was  proposed  by  Coons  in  MAC-TR-41.  The  particular  implementation 
discussed  and  the  derivation  of  the  necessary  algebra  and  matrices  to 
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make  it  work  are  the  author's  responsibility,  as  is  the  error  analysis. 

The  subsection  on  surface  construction  utilizing  interior  curves  was 
suggested  to  the  author  in  considerable  detail  by  Coons.  To  the  best 
of  my  knowledge,  the  rest  of  the  material  is  original.  Most  subsections 
of  Sections  II  and  III  were  presented  in  almost  their  present  form  as 
two  separate  papers  to  the  1969  Spring  Joint  Computer  Conference  [5,  13]. 

The  remainder  of  this  introductory  section  presents  the  context 
in  which  this  work  has  arisen  and  furnishes  a  uniform  notation  and 
format  for  the  homogeneous  coordinate  geometry.  The  subsection  on 
notation  should  in  particular  be  held  firmly  in  mind  when  reading  the 
mathematics  of  the  surface  formulation. 

1.2  ANTECEDENTS  AND  HISTORICAL  RELATIONS 

The  existence  and  direction  of  this  research  have  been  influenced 
by  many  factors,  both  intentional  and  accidental.  The  organization  of 
this  thesis  is  approximately  in  chronological  sequence,  but  only  because 
the  topical  dependencies  guided  the  research  in  the  same  direction  as 
seems  best  for  presenting  the  material.  In  this  brief  subsection  we  con¬ 
sider  the  place  of  this  work  in  its  historical  context. 

In  the  initial  stages  the  research  was  concerned  purely  with  the 
problems  of  using  the  proposed  Harvard  Three-Dimensional  Display  [21] 
to  draw  curves  as  well  as  straight-line  figures.  Taking  our  direction 
from  the  earlier  work  of  Blatt  and  Roberts  at  the  M.  I.  T.  Lincoln 
Laboratory  [4,  16]  using  rational  quadratics  it  was  suggested  that  since 
the  equipment  was  to  use  a  4  X  4  matrix  we  ought  to  consider  adding  a 
cubic  term  for  additional  flexibility.  The  results  presented  in  Section  II 
on  curves  indicated  that  this  was  indeed  a  promising  approach  and  led  to 
the  inclusion  of  appropriate  logic  in  the  hardware  to  implement  the  for¬ 
ward  difference  scheme  in  the  display. 

After  having  in  the  course  of  investigating  the  properties  of  the 
curve  formulation  come  to  some  basic  understanding  of  the  simplicity  of 
the  parametric  polynomial  approach,  we  suggested  expanding  the  poly¬ 
nomial  to  two  variables,  thereby  generating  surfaces.  Having  demonstrated 
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by  construction  the  generation  of  the  sphere  and  the  cylinder,  we 
abandoned  immediate  research  on  curves  to  proceed  more  directly 
with  surfaces.  Some  time  was  spent  in  becoming  familiar  with  Coons' 
earlier  work  [6,  7],  but  it  was  decided  that  as  far  as  this  research  was 
concerned,  adapting  his  more  general  formulation  to  homogeneous 
coordinates  would  lose  the  efficiency  and  simplicity  of  this  approach. 

The  notation  has  however  been  gratefully  borrowed. 

It  soon  became  clear  that  the  research  would  be  progressing 
more  rapidly  than  the  construction  of  the  three-dimensional  display. 

We  thus  changed  the  character  of  the  research  slightly  to  make  it 
more  mathematical  rather  than  directed  to  the  utilization  of  this 
specific  equipment.  This  change  is  what  has  limited  experimental 
implementation  to  the  rudimentary  demonstration  programs  of 
Subsections  4.1  and  4.2  rather  than  to  a  more  complete  system,  such 
as  proposed  in  Appendix  II.  If  we  had  proceeded  in  that  direction, 
moreover,  the  resulting  effort  would  probably  not  have  been  as  general 
as  is  now  proposed. 

The  evidence  indicates  that  other  work  in  this  field  --  either 
derived  from  Coons'  formulation  or  from  various  implementations  of 
interpolation  in  three  dimensions  --as  performed  by  the  airframe  or 
automobile  industries  is  not  strictly  relevant  to  this  formulation.  Even 
though  we  have  lost  immediate  contact  with  the  three-dimensional  hard¬ 
ware,  we  are  concerned  primarily  with  an  efficient,  flexible,  relatively 
simple  way  of  drawing  three-dimensional  curved  surfaces  subject  to 
certain  constraints.  This  is  to  enable  an  interactive  approach  to  the 
specification  of  surfaces  without  requiring  excessive  computation  power. 
The  existing  interactive  systems  of  which  we  are  aware  are  only 
approximately  so,  limited  by  ad  hoc  procedures  relevant  to  the  particu¬ 
lar  problem  area  (e.g.,  airframe  specification  and  analysis)  rather  than 
to  a  more  general  class  of  problems  encompassing  a  wide  variety  of 
surfaces.  Strictly  speaking,  these  results  are  not  at  all  relevant  to  the 
problem  of,  for  instance,  fitting  a  smooth  surface  through  a  set  of 
measured  points  on  a  clay  model  in  order  to  drive  a  machine  tool.  These 
surfaces  could  be  used  for  this  kind  of  three-dimensional  interpolation, 
but  only  rather  clumsily. 
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1.3  NOTATION 


We  will  be  talking  about  surfaces,  represented  as  tensors,  curves, 
represented  as  matrices,  or  points,  represented  as  either  three- 
dimensional  (ordinary  coordinates)  or  four-dimensional  (homogeneous 
coordinates)  vectors.  A  vector,  usually  a  point  in  homogeneous  coordi¬ 
nates,  will  always  be  denoted  by  a  vector  arrow,  for  example,  V. 

Where  relevant,  a  four-dimensional  vector  will  be  represented  by  an 
upper-case  letter  and  a  three-dimensional  vector  by  a  lower-case  letter. 
The  array  representing  a  point  or  a  curve  may  be  independent  or  it  may 
be  a  portion  of  a  higher  dimensional  array  representing  a  higher  order 
object. 


Subscripts  will  be  used  to  denote  either  components  of  the  array 

(tensor,  matrix,  or  vector)  or  to  indicate  partial  derivatives  with 

respect  to  the  parameters.  Components  of  vectors  will  not  have  a 

vector  arrow,  although  vector  components  of  a  matrix  will,  for  example, 

the  vector  of  the  matrix  A^.  Integer  subscripts  .,  ^  or  explicit 

numeric  subscripts  will  be  used  for  components  in  some  cases  while  the 

symbols  ^x,  ^y*  ^z,  h  be  usec^  when  it  is  desirable  to  emphasize  the 

spatial  coordinates.  The  subscripts  ,  will  naturally  refer  to  the 

3  3  u  v 

partial  derivatives  7^  ,  —  .  The  components  of  a  point  in  three 

dimensions  are  indicated  by  subscripts  ,  ,  --  for  example,  p  . 

x  y  z  x 

This  implies  either  true  three-dimensional  data  or  a  division  to  remove 
the  homogeneous  scale  factor. 


When  a  tensor  describing  a  surface  is  used  without  any  subscripts 
for  components  it  is  to  be  treated  not  as  an  array  of  scalars  but  rather 
as  a  matrix  whose  elements  are  vectors.  Multiplication  of  such  a 
matrix  of  vectors  by  a  matrix  of  scalars  is  to  be  interpreted  in  the 
standard  way,  with  the  summation  being  a  vector  sum  and  the  individual 
multiplications  being  the  product  of  a  scalar  and  a  vector. 

When  convenient,  the  standard  convention  of  summing  over  repeated 
indices  will  be  used,  with  the  proviso  that  indices  appearing  on  both  sides 
of  an  equation  will  not  be  summed  but  indicate  a  running  index.  In  some 
cases  an  indicial  quantity  will  be  used  as  an  actual  exponent  --  context  will 
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make  it  clear  when  that  occurs;  in  any  case,  if  such  an  exponent  is 
repeated  in  an  equation  it  is  to  be  treated  according  to  the  summation 
rule.  Such  summation  and  running  indices  will  take  on  the  values 
0,  1,  2,  3  for  virtually  all  of  the  thesis. 

A  curve  ~!§(u,  v),  u=a,  a=constant,  1,  on  a  surface  ^(u,  v) 

will  be  denoted  by  a  or  by  J?av  for  brevity.  A  point,  ^(u,  v),  u=a, 
v=b,  on  the  surface  will  be  denoted  by  u  a>  v  °  or  by  ^a^.  In  some 
cases  these  symbols  will  be  further  abbreviated  by  dropping  the 
symbol  S  for  the  surface  and  writing  just  the  values  of  the  parameters  -- 
see  subsection  3.3.  In  most  cases  it  does  not  matter  whether  such 
notation  is  taken  to  represent  the  object  itself  or  the  array  describing 
the  object;  context  should  make  the  situation  clear.  The  geometric  and 
mathematical  bases  in  which  such  an  array  represents  an  object  will 
either  be  irrelevant  or  clearly  specified. 

1.4  HOMOGENEOUS  COORDINATES 

Taking  the  lead  from  the  19^^  century  projective  geometrists  [1], 
computer  display  programmers  have  recently  adopted  the  use  of  a 
homogeneous  coordinate  representation  of  geometrical  data.  This 
formulation  is  especially  relevant  when  talking  about  perspective 
images.  As  is  to  be  expected,  most  of  the  results  derived  by  the 
mathematicians  do  not  seem  relevant  to  computer  applications  for  they 
normally  involve  deductive  rather  than  constructive  proofs.  This  sub¬ 
section  summarizes  a  few  of  the  definitions,  conventions  and  formulae 
in  homogeneous  analytic  projective  geometry  needed  for  an  understanding 
of  the  surface  formulation.  Most  of  this  material  has  been  printed  else¬ 
where,  but  is  not  readily  available  [14,  15,  17,  24], 

Point.  A  point  (in  three  dimensions)  is  represented  by  a  non-zero 
vector  of  four  components.  Two  points  ?,  C§  are  to  be  regarded  as  the 
same  point  if  and  only  if  they  are  linearly  dependent;  that  is,  if  and  only 
if  there  exists  a  non-zero  constant  a  such  that  ?  =  uC$.  The  ordinary 
three-dimensional  point  [x,y,  z]  can  and  will  be  represented  in  the 
homogeneous  form  as  [x,  y,  z,  1];  hence  any  point  [a,  b,  c,d],  d^O,  is 
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cL  b  c 

equivalent  to  the  three-dimensional  point  [-j  ,  -j-  ,  J .  This  represen¬ 
tation  will  be  indicated  by  the  notation  [hx,  hy,  hz,  h]  =  hv  =  h[x,  y,  z,  1]  , 
v  =  [x,  y,  z,  1]  for  the  homogeneous  coordinates  of  a  point  V,  where  by 
implication  we  take  the  quotients  ^  ^  to  find  the  three- 

dimensional  coordinates.  The  assumption  is  that  once  a  point 
[hx,  hy,  hz,  h]  has  been  obtained,  something  --  hardware  or  software  -- 
will  perform  the  necessary  division  by  the  homogeneous  coordinate.  In 
particular,  whenever  we  talk  about  a  point  to  be  displayed,  this  division 
must  be  performed. 

Line:  Three  points  ?,  ($,  R  are  collinear  if  and  only  if  they  are 
linearly  dependent.  The  set  of  points  ?  on  the  line  through  two  points 
Q,  R  can  thus  be  generated  parametrically  by  P  =  aQ  +  (l-a)R,  or,  in 
full  generality,  by  ?  =  +  /3lt,  a,  not  both  zero.  Two  points  on  this 

line,  and  ~  +  ^2^  are  equivalent  if  and  only  if 

the  vectors  and  [ ]  are  linearly  dependent,  that  is,  in  this 

case,  proportional. 

Plane:  Four  points  P,  Q,  R,  S,  are  coplanar  if  and  only  if  they 
are  linearly  dependent.  Hence  the  plane  through  three  points  P,  <5,  ^ 
is  spanned  by  oP  +  /3C§  +  7K,  a,  fi,  7  n°t  all  zero.  Furthermore,  for 
the  dependence  +  /3C$  +  7^  +  6!?  =  0  to  hold,  we  must  have 
I  P  Q  R  S|  =  0,  or,  expanding  by  minors  on  the  fourth  column,  S  •  T  =  0 
(vector  dot  product)  where 


P2  ^2  R2 

piQi  Ri 

piQiRi 

piQiRi 

P3  Q3  R3 

3 

P 3  ^3  R3 

3 

P2  ^2  R2 

3 

P2  ^2  R2 

P4  Q4  R4 

P4  Q4  R4 

P4  Q4  R4 

P 3  ^3  R3 

The  relation  is  •  'T'  =  0  is  thus  the  equation  of  a  plane  and  the 
vector  T  represents  the  plane.  Among  many  results  it  can  be  shown 
that  two  planes  T  and  U  are  equivalent  if  and  only  if  the  vectors  T 
and  U  are  proportional. 


11 


1.5  TRANSFORMATIONS 


In  ordinary  three  space,  a  non-singular  matrix  transformation  of 
the  form  =  P^T.^  (T  a  3  X  3  matrix)  is  an  affine  transformation,  pro¬ 
ducing  only  rotation  and  scaling.  In  the  special  three-dimensional  space 
of  homogeneous  coordinates,  such  a  transformation  (where  T  is  now 
4X4)  is  called  a  projective  transformation  which  in  addition  to  rotation 
and  scaling  performs  translation  and  a  perspective  transformation, 
hence  its  applicability  to  computer  graphics.  This  transformation  is 
derived  as  follows  (see  Figure  1.5.1  and  Reference  [22]): 

Let  an  observer  be  at  the  origin  5  of  a  coordinate  system  and  let 
a  display  screen  S  of  size  2  scope  units  by  2  scope  units  be  at  position 
z  =  cot  a,  where  a  is  half  the  angle  subtended  by  the  screen.  The  x 
coordinate  xg  of  the  intersection  with  the  screen  of  a  ray  from  a  point 
P  to  the  observer  --  perspective  projection  from  O  of  P  onto  S  --  is 


Similarly, 

ys  '  HS  cot  “ 


and  let  us  define  an  inverse  distance  as 
h 

Zs  hz 


Then  the  point  represented  by 


P 

s 


hx  , 
r—  cot  a, 
hz  * 


hz 

hz 


cot  a, 


_h 
hz  ’ 


1 


gives  as  its  x  and  y  coordinates  the  location  on  the  screen  at  which  to 
draw  the  perspective  projection  of  the  point  P  and  as  its  z  coordinate 
a  value  monotonically  related  to  distance,  suitable  for  intensity  modulation. 
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Fig.  1.5.1  Perspective  Projection 
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P  is  equivalent  to  (hz)P  =  [hx  cot  cr,  hy  cot  o',  h,  hz]  = 
S  S 


=  [hx,  hy,  hz,  h] 


c  ot  a  0  0  0 

0  c  ot  a  0  0 

0  0  0  1 

_  0  0  1  0 


and  we  have  a  matrix  which  performs  a  perspective  projection  on  a 
point  represented  in  homogeneous  coordinates. 

Translation  is  achieved  by  a  matrix  of  the  form 

”l  0  0  0~ 

0  10  0 

0  0  10 

xt  yt  zt  1 


and  rotation  by  a  matrix  of  the  form 


0  " 

R 

0 

0 

0  0  0 

1 J 

where  R  is  an  ordinary  three-dimensional  rotation  matrix.  Typically, 
several  such  matrices  --  rotation,  translation,  and  projection  --  will 
be  applied  in  sequence  to  arrive  at  a  single  compound  transformation 
expressed  as  the  matrix  product  of  the  separate  matrices. 

Observe  that  given  any  matrix,  M,  which  performs  some  trans¬ 
formation,  a  matrix  cyM,  a  4  0,  performs  the  same  transformation, 
since  cKP  X  M)  =  P  X  cyM  implies  P  X  M  and  P  X  oM  are  proportional, 
for  any  point  P.  In  particular,  the  identity  transformation  becomes  al , 
a  4  0. 
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1.6  DISPLAY  OF  LINES 


We  have  decided  to  use  homogeneous  coordinates  for  the  repre¬ 
sentation  of  points  in  three-dimensional  space  with  an  intent  to  display 
such  points  in  perspective  on  a  two-dimensional  screen.  In  general  we 
will  be  talking  about  drawing  straight  lines  rather  than  merely  isolated 
points.  Each  such  straight  line  will  be  represented  by  a  starting  point 
P  and  a  finishing  point  (^,  both  ?  and  ^  represented  in  homogeneous 
coordinates.  To  display  a  perspective  view  of  a  line  we  apply  the  per¬ 
spective  transformation  to  the  end  points  of  the  line  and  display  a  two- 
dimensional  vector  between  the  projections  of  the  endpoints. 


Let  R  =  PT  and  S  =  QT  where  T  is  the  desired  perspective 


transformation  matrix.  Then 


Rx  R, 
FT.  ’  R 


4  J 


and 


4  J 


are  the  pro¬ 


jected  endpoints  of  the  line  on  a  display  screen  whose  range  is  -1<x<1, 
R  3  S3 

- 1  <  y  <  1  •  ■5-  and  -s—  are  for  the  startpoint  and  endpoint,  respectively, 

R4  b4 

the  reciprocal  of  the  distance  from  the  point  to  the  observer,  .  If  we 
apply  this  third  coordinate  to  the  intensification  control  on  a  display,  we 
will  obtain  a  picture  in  which  the  brightness  varies  approximately 
inversely  with  distance,  depending  on  the  response  of  the  display  to  the 
intensification  control.  (The  Adage  commercial  display  uses  a  similar 
technique  for  depth  indication.)  As  the  vector  generator  drawing  the 
line  from  ?  to  ^  moves  the  beam  across  the  scope,  we  would  want  it 
to  change  the  intensity  in  a  correlated  manner  from  the  value  at  the 
start  to  the  value  at  the  finish.  As  -*•  «>,  (z  -►  0),  we  would  have  to 
limit  the  brightness  to  some  fixed  level;  correspondingly  we  might  not 
want  =  0  (z  =  °°)  to  be  represented  as  absolutely  no  intensity,  but  as 
some  just  barely  visible  minimum. 


Notice  the  following  two  difficulties  in  a  naive  application  of  the 
above: 


Rf  Rg 

1.  What  if  one  or  more  of  the  ratios  ■=-  ,  =- 

^4  ^4 


■q— ■  is  greater 
b4 


than  one  in  absolute  value? 
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R3  S3 

2.  What  if  one  or  both  of  the  intensities  ,  o—  is  less  than  zero? 

K4  b4 

Condition  1  means  that  the  point  of  the  line  causing  the  difficulty  is 
outside  the  viewing  window  defined  by  the  edges  of  the  display. 

Condition  2  means  that  the  appropriate  point  of  the  line  is  behind  the 
observer  (z  <  0)  and  should  not  be  visible. 


The  Clipping  Divider  [19]  associated  with  the  Harvard  Three- 
Dimensional  Display  [21]  was  designed  to  solve  the  above  problems  by 
determining  digitally  before  the  division  whether  the  endpoints  of  the 
line  lie  outside  the  viewing  window,  and,  if  they  do,  to  compute  the 
intersection  of  the  line  with  the  window  in  order  that  a  vector  may  be 
drawn  between  two  points  within  the  display  area.  (See  Figure  1.6.1.) 
The  Clipping  Divider  does  not  compute  the  term  and  does  not  deal 
with  the  problem  of  intensity  modulation.  It  assumes  that  the  homo¬ 
geneous  coordinate  h  of  the  original  points  is  positive;  it  thus  only 
uses  the  signs  of  R^  or  to  tell  whether  the  point  is  in  front  of  or 
behind  the  observer  and  does  not  care  how  far.  It  does  compute  the 
intersection  for  a  line  which  passes  from  in  front  of  the  observer  to 
behind  him,  or  vice  versa. 

In  the  software  implementation  of  this  clipping  or  windowing 
process,  a  modified  version  of  a  routine  written  by  D.  Cohen,  we  also 
do  not  provide  for  intensity  modulation  and  avoid  the  division.  To 
include  these  would  have  slowed  the  program  too  much,  and,  would  not 
have  given  a  good  result  for  the  effort  since  the  intensity  control  on  the 
DEC  type  340  display  used  is  not  suitable  for  the  fine  variation  in 
intensity  needed.  Conversations  with  A.  R.  Forrest,  Joint  Computer- 
Aided  Design  Group,  University  of  Cambridge,  indicate  that  a  half 
dozen  or  so  intensity  levels,  appropriately  spaced,  selected  by  a  soft¬ 
ware  mapping  from  t^-  gives  a  realistic  appearance;  even  this  level  of 
intensity  control  was  not  readily  possible  with  the  equipment. 
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Display  Screen 


Entire  Line 
In  three  dimensions) 


FIG.  t.6.1  CLIPPING  OF  A  LINE 
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SECTION  II 


CURVES 


2.1  INTRODUCTION 

Now  that  the  computer  displays  with  vector  and  character 
drawing  capabilities  are  becoming  a  common  form  of  on-line  graphic 
output  device,  interest  has  turned  toward  providing  a  curvilinear 
display  capability  [4,  10,  16,  20].  Naturally,  curves  can  be  drawn  as  a 
series  of  short  vectors  stored  in  display  memory,  but  the  goal  of 
current  research  is  to  provide  hardware  which  draws  a  curve  by 
generating  such  vectors,  or  perhaps,  continuous  beam  motion,  from 
a  concise  specification  of  a  curve  segment. 

In  this  chapter  we  consider  the  properties  of  a  particular  class 
of  curves  suitable  for  such  a  hardware  curve  generator.  (See 
Figures  2.1.1  and  2.1.2  for  examples.)  As  it  happens,  the  algorithms 
for  generating  successive  points  on  the  curve  can  be  implemented 
efficiently  in  software. 

Mathematically  we  present  this  material  in  as  general  a  form  as 
seems  appropriate  for  the  particular  topic  under  discussion.  Practi¬ 
cally  we  think  primarily  of  two-dimensional  curves  or  three-dimensional 
curves  projected  into  two  dimensions. 

We  will  consider  the  family  of  curve  segments 

<£  =  {V  =  V(t)  ,  0  <  t  5  1} 

where  each  component  of  V  is  a  polynomial  in  t.  We  restrict  the  dis¬ 
cussion  to  polynomials  only  because  of  the  simplicity  of  polynomial 

N 

evaluation.  We  initially  consider  these  curves  in  R  ,  the  space  of  real 

N 

vectors  of  dimension  N;  however,  our  interest  lies  in  P  ,  the  per- 

N 

spective  space  of  dimension  N-l,  generated  from  R  by  dividing  each 
component  of  V  by  the  last  component. 

Let  <j>  be  the  family  of  curve  segments  defined  by  polynomials  of 
degree  not  exceeding  m.  The  higher  m  is,  the  more  general  is  the 
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Fig.  2.1.1  Example  of 


I  X 

a  Curve 
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x 


Fig.  2.1.2  A  Closed  Curve 
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family  <f>  ,  but  generating  each  point  on  any  curve  is  more  complex. 

The  greater  generality  of  for  larger  m  enables  the  curve  to  satisfy 
more  conditions,  but  correspondingly  requires  more  conditions  to 
uniquely  define  the  curve.  $ 2  includes  straight  lines  and  conic  sections 
only.  Therefore  if  V  e  (f> V  cannot  cut  itself  and  cannot  have  inflection 
points.  (which  includes  as  a  proper  subset)  contains  curve 
segments  which  may  have  inflection  points  and  which  may  cut  themselves, 
as  shown  in  Figures  2.1.3  and  2.1.4. 

2.2  ITERATIVE  GENERATION  OF  CURVE 

In  order  to  display  a  curve  we  wish  to  compute  successive  points 
along  the  curve,  forming  a  straight-line  approximation  to  the  curve  by 
connecting  adjacent  points  with  short  display  vectors.  Let 

V(t)  =  [  Vj(t),  v2(t),  v3(t),  V4(t)  ] 

be  the  equation  of  the  curve  in  homogeneous  coordinates.  Each  of  the 
components  V.(t)  of  V(t)  is  a  polynomial  in  t.  Let  P(t)  represent  any  of 
these  polynomials.  Although  the  successive  points  { t^,  t2,  t^,  ...,  t-^} 
could  be  any  such  that  0  S  tg  <  t.  <  ...  <  ^  <  t^  S  1,  we  will  choose  to 

space  the  points  equally  along  the  parameter  Define  Pn  =  P(t  )  =  P(n6) 
where  6  =  ^7.  We  could  compute  the  points  by  evaluating  P  for  each  t  , 
but  such  an  approach  would  be  highly  inefficient 

We  choose  therefore  to  use  a  finite  difference  method.  Let  e  be  a 
scale  factor  (to  be  discussed  shortly).  Let  f  be  any  function.  We  define 
two  difference  operators  with  respect  to  the  difference  6‘ 

the  forward  difference:  Af(x)  = -^  [  f(x+ 6)  -  f(x)] 

the  backward  difference:  Vf(x)  = -^  [  f(x)  -  f(x  -  6)] 

Am  and  Vm  are  defined  recursively  where  A^f  =  V^f  =  f .  Note  that 

1  m  t 

—  [P  ...  -P  1  =  AP  =  VP  ,  1  and  that  if  P(t)  =  7  a,  t  then 
eLn+lnJ  n  n+1  wz_/k 

k=0 

AmP  =  VmP  =  —  m!  a  for  any  value  of  t.  To  be  specific,  we  will 

e 
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Fig.  2.1.3  Curve  with 
Inflection  Point 


# 

Fig.  2.1.4  Curve  which  Crosses  Itself 
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choose  to  let  m=  3,  that  is,  we  will  talk  about  cubic  polynomials. 
If  we  define  P  =  [x3x2xl]A,  A  e  r\  then 


and 
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It  is  easy  to  see  that 
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A2P 

a3p 


n+1 


n+1 


n+1 


which  we  will  abbreviate  to  F 
backward  differences. 
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=  SF  or  as  F  =  SnF„  and  for  the 
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abbreviated  as  B 


=  TB  or  B  =  TnBn  . 
n+1  n  n  0 
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where  we  identify  as  the  operation  of  scaling  A1+*P  by  e  and  adding 
the  scaled  value  to  A1?  (or,  V1+^P  to  VXP).  These  scaling  and  addition 
operations  can  be  performed  serially  in  the  order  A^,  A^,  A2  for  forward 
differences  and  in  the  opposite  order  A2,  A^,  A^  for  backward. 

In  both  schemes  there  are  only  three  additions;  however,  in  the 
forward  difference  scheme,  each  new  value  A*P  .  depends  on  old 


1  i+ 1  n+ 1 

values  AP  ,  A  iP  whereas  for  the  backward  difference  scheme  each 

ni  n  3 

new  value  V  Pn+^  except  for  V  Pn+^  (which  is  constant)  depends  on 

both  the  old  value  V1?  and  the  new  value  of  V1+*P  .  , . 

n  n+1 


Thus  all  the 

operations  of  the  forward  difference  scheme  could  be  done  in  parallel, 

while  those  for  the  backward  must  be  done  serially.  If  the  forward 

difference  scheme  is  done  in  parallel,  the  old  values  A1?^  cannot  be 

destroyed  until  all  the  new  values  are  obtained;  if  either  scheme  is 

done  serially,  each  A*P  (  or  VXP  )  can  be  replaced  by  A^P  .  . 

(or  V  P  as  soon  as  it  is  computed.  Using  serial  programming,  the 

backward  difference  scheme  is  slightly  faster  since  the  result  V 

of  one  operation  is  immediately  available  to  compute  the  next  value 

V*  Ip^.  and  need  not  be  fetched  from  memory;  in  the  forward  scheme 

both  the  AXP  .  and  A1+^P  .  needed  to  compute  A*F  must  be  fetched. 

n-1  n-1  r  n 
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In  all  the  above  we  replace  the  polynomial  P  by  the  vector  V  and 
we  have  an  efficient  iterative  method  for  generating  the  successive 
points  V(n6)  on  the  curve,  keeping  the  4X4  matrix  of  differences 
(A1V(n6)  or  V1V(n6))  in  memory  or  fast  registers. 

If  we  represent  V  =  [t  t  1 1]  A,  Aa4X4  matrix,  then  the  initial 
set  of  differences  is  given  by 
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or  by 
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The  scale  factor  e  is  to  be  chosen  to  minimize  the  accumulated 
rounding  error  in  computing  V(t  ).  We  assume  all  numbers  are  stored 
as  fixed  point  fractions.  Since  we  are  dealing  in  homogeneous  coordi¬ 
nates,  we  can  multiply 

v(or 

A  V(0) 

A2  V(0) 

.  A3V(0)  _ 

by  a  normalization  constant  a  in  order  to  prevent  any  of  the  components 
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of  A1:V(t.)  from  exceeding  one  in  magnitude.  The  normalization  a 

3 

should  be  chosen  as  large  as  possible  in  order  that  as  much  precision 
as  possible  be  used.  We  will  analyze  the  forward  difference  algorithm 
in  terms  of  a  single  polynomial  P  rather  than  the  vector  V,  for 
simplicity. 


We  can  write  the  forward  difference  scheme  in  matrix  form  as 
Fn  =  SnFp  =  (I+eH)nFp  where  in  general  (H^) _  =  or  for  k  =  1,  2,  3, 
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0 
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and 


Expanding  (I+eH)n  =  ^  e^H*,  we  get 

i 

l  (i)^% 

i=0 


For  the  backward  difference  scheme,  we  have  similarly 
§n  =  (1+  eH  +  e3H3  +  e3H3)n  Now  we  can  write  (1+  eH  +  e3H3  +  e3H3)  = 

=  (I  -  e^H^)(I  -  eH)  *  =  (I-eH)  *  since  =  0.  Expanding  formally,  this 

3 

gives  (1+ eH  +  e3H3  +  e3H3)n  =  (I  -  eH)  n  =  ^  where  the 

i=0 

4 

summation  stops  at  i  =  3  since  H  =0. 

By  a  standard  identity  [12,  p.  56],  ^  =  (-1)1  (n+|  ^  and  we 

have  for  the  backwards  differences. 


VJP 


n 


l  (n+['1)e1viflP0 

i=0 
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We  will  analyze  the  forward  difference  scheme  only,  in  order  to 
determine  optimal  values  for  e  and  6.  The  analysis  of  the  backwards 
scheme  would  be  similar,  but  with  additional  algebraic  complexity. 


Define  A^P  -  A^P  =  e^  where  A%  is  the  correct  value,  A^P  the 
n  n  n  n  J  n 

computed  value,  and  e^  the  accumulated  error.  At  any  step  of  the 


n 


process,  we  will  have  A^P  =  p^  +  A^P  .  +  eA^4  ^P  .  where  p^  is 
r  >  n  Kn  n-1  n-1  Kn 

the  roundoff  error  introduced  in  the  computation.  Define  e^  =  p^  ,  that 

is,  the  initial  error  in  A^P^  is  due  to  roundoff. 

Subtracting  the  original  difference  equation  for  A^Pn,  we  obtain 

AJp  -  a^P  =  pJ  +  (A^P  ,  -AJP  1 )  +  e(A^+1P  ,  -AJ>1P  J  or 

n  n  Kn  n-1  n-1  n-1  n-1 

=  pJ  +  1  +  eeJ+1.  . 

n  Kn  n-1  n-1 
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L  n 

-  n 

we  have  in  vector  form 


En  =  (I+eH)Sn_1  +  Rn  =  (I+eH)[(I+eH)En_2  +  +  Rn  = 

=  (I+eH)2S  .  +  (I+eH)^  .  +  R  = 

n-2  n-1  n 

n-1  .  n 

=  ]]  (I+eH)1^n_i  +  (I+eH)nEQ  =  V  (i+eH)1^ 

i=0  i=0 

n  3  . 

panding  (I+eH)1  as  before,  we  have  =  Y  ^  .  or 


Ex 


4 - 1  X  (ic)-k[s  4+k4-J  - 1  l  (i) «k4S  • 

i=0  k=0  L  6.  J  i=0  k=0 
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Now  let  us  assume  that  the  rounding  errors  are  given  by  inde- 

2 

pendent  random  variables  with  the  same  mean  and  variance  a  .  We 

q  q  v  O  v 

will  also  make  use  of  the  fact  that  A^P  =AF  .  (p  =0  for  n  >  0), 

n  g  n-1  Kn  ’ 

since  no  computation  has  to  be  done  on  the  A  term.  In  other  words. 


and 


Thus  E3 

i=0  k=0 


= 


[i  1! 


By  a  standard  identity  [12,  p.  54],  J  .  we  get 


i=0 


»el  =  lu 0 


r2f 

k=0  J 

2  j 

v  n J 


To  compute  the  variance  a  e^  recall  that  for  any  sum  of  inde¬ 
pendent  random  variables, 

2  ~  v  2  2 


=  )  a.x.  ,  a  x  =  )  a.  a  x. 

Li  ii’  Li  i  i 


so  that 


2  j  2 

cr  e  =  cr„ 
n  0 


-lo  t2<3'j)- 


By  a  not-so-standard  identity  [11,  p.  106], 
n  .  9  k 

l  ' 

i=0  '  j?=0 


.2  (L)  ■  l 
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2-j 


k 

z 

£=0 


.  2  j  2  v  2k  Y  /  ^i/n+i+lWn+n  ,  /  n  \  2(3-j) 

ami  a  en  “  CTo  |_  Z  e  L  ^  U+i+1  )  U-  J  +  U-j  J  E  J 

k  0  ^ 


The  expanded  expressions  for  the  combinatorial  sums  are  as 
follows: 
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Notice  that  all  the  coefficients  of  the  powers  of  n  are  of  the  same 
magnitude  in  each  expression;  we  can  thus  reasonably  approximate 
these  error  estimates  by  replacing  each  such  expression  by  its  first 
term,  the  approximation  being  asymptotically  valid  for  n  »  1.  The 
resulting  simplified  expressions  for  the  errors  are  thus 
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For  ne  -»■  0,  these  expressions  take  on  the  limiting  values  /uen  =  /u _n, 
o“en  =  62n  or  ae^  = 
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The  interpretation  of  these  results  is  that  most  of  the  roundoff 

error  occurs  as  introduced  in  the  last  stage  (V  +  sAV  -►  V  .  ,).  A 

second  conclusion  is  that  it  is  significantly  worth  it  to  use  an  unbiased 

rounding  scheme  --  such  as  rounding  on  the  basis  of  the  first  bit  lost 

in  the  multiplication  by  e  --  rather  than  mere  truncation.  This  is 

because  the  mean  of  the  rounding  errors  propogates  much  faster  than 

the  standard  deviation.  In  fact,  if  truncation  errors  permit  only  on  the 

order  of  m  steps  to  be  taken,  correct  rounding  will  permit  on  the 
2 

order  of  m  steps. 

If  the  roundoff  errors  (Mq,  Oq)  in  the  initial  A^Pq's  are  larger 
than  the  errors  occurring  in  the  iterations  (Mq,  ct^)  --  as  would 
normally  be  the  case  since  the  A^P^'s  are  the  results  of  several 
computations  --  the  terms  for 


k  l+k 
e  PJQ 


will  have  to  be  specifically  added.  This  would  not  change  the  error 

k  k 

significantly,  since  it  would  only  add  terms  of  the  form  n  e 

whereas  the  corresponding  terms  already  are  of  the  form 
k  k 

n(n  e  a  factor  of  n  larger,  and  similarly  for  cr. 

The  absolute  error  derived  above  is  independent  of  the  actual 

polynomial  computed.  The  relative  error,  which  measures  the  real 

0 

e 

accuracy  of  the  computation,  must  be  computed  as  —  and  is 

n 

dependent  on  P(x).  If  we  pick  the  normalization  factor  a  as 

0 


«  =  1  /( 


maxi 


AJP, 


!)• 


the  relative  error  is  thus 


n 


aP 


Figure  2.2.1  is 


n 


ue 


.0 

n 


a  graph  of  the  normalized  relative  mean  error  defined  as  — -  for 

6  ^  aP  unn 

n  0 

3  2  2 

the  polynomials  x  +x  +x+l,  x  +x+l  and  x+ 1,  three  typical  cases. 
These  graphs  are  for  n=64.  The  graphs  for  other  n  and  for  the  normal¬ 
ized  standard  deviation  of  the  error  would  be  similar,  both  in  shape  and 
value. 
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P=  x3+x2+x  +  i 


3.  2-2.1  NORMALIZED  RELATIVE  MEAN  ERROR 
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3  2 

For  the  case  P=x  +x  +  x+  1,  the  maximum  occurs  at  x=  1, 
namely. 


f  p  1 
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For  (ne)  <  1,  the  dominating  term  is  - ^  while  for  (ne)  >  1,  the 

(ner 

maximum  becomes  a  constant  4.  The  relative  errors  are  then  given  by 

dividing  this  maximum  by  the  value  assumed  at  P  and  multiplying  by 

n  0  0 

the  absolute  error.  Using  the  asymptotic  values  for  ue^  and  cren,  we 
have 
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The  proper  compromise  is  to  choose  ne  =  1,  that  is,  n  =  —  and  6 
as  confirmed  by  the  exact  data  in  Figure  2.2.1 


=  e, 


The  final  conclusions  on  selecting  n  and  e  are  thus: 


1.  Conclusion:  Select  n  as  small  as  possible  consistent  with  the 
number  of  line  segments  needed  to  accurately  represent  the  curve. 

Reason:  The  mean  of  the  error  increases  linearly  with  n  and 
the  standard  deviation  as  the  square  root  of  n. 

2.  Conclusion:  Select  e  =  —  and  certainly  not  outside  of  the  range 
12  n 

2^  <  e  <  — .  e  may  be  chosen  slightly  less  than  optimal  if  the  curve  is 
relatively  smooth. 

Reason:  The  errors  vary  roughly  with  the  square  of  (ne)  or 
depending  on  which  is  worse,  and,  for  (ne)  <  1  are  worse  the  higher  the 
degree  of  the  polynomial. 
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Notice  that  the  apparent  error  in  generating  higher  order  curves 
is  particularly  bad  --  n  must  be  large  in  order  that  the  line  segment 
approximation  be  good,  and  the  presence  of  cubic  terms  heightens  the 
error  caused  by  scaling. 

In  applying  this  analysis  to  actual  curve  generation,  we  must 
recall  that  the  scaling  is  applied  to  all  components  of  V(t)  equally; 
hence  a  cubic  term  in  one  increases  the  relative  error  in  all. 

Of  course  these  are  statistical  estimates.  However,  if  /Uq  is  taken 
as  the  absolute  value  of  the  maximum  rounding  error,  then  qPn  will  be 
an  upper  bound  on  the  magnitude  of  the  absolute  accumulated  error.  On 
the  other  hand,  assuming  that  the  rounding  errors  are  distributed 
randomly  may  not  be  valid  --  certain  polynomials  with  coefficients 
expressed  exactly  with  few  enough  bits  and  drawn  with  few  enough 
segments  can  be  computed  exactly.  Experiments,  however,  do  indicate 
that  the  relation  of  rounding  error  to  (ne)  is  very  strong,  as  predicted 
by  theory. 

These  results  were  for  generating  curves.  For  drawing  surfaces  -- 

using  a  double  iteration  --  the  absolute  error  is  approximately  twice  as 

large  and  the  scaling  the  square  of  the  scaling  in  the  case  of  curves. 

3  3 

Thus  for  P(x,  y)  =  x  y  ,  the  relative  error  of  the  mean  is  on  the  order  of 
2  •  36  •  n,  or,  for  n=32,  about  2^,  also  verified  by  experiment. 

A  possible  variation  in  this  scheme  would  be  to  define 

Af(x)  =  —  (f(x+6)  -f(x)) 
el 

A^f(x)  =  —  (Af(x+6)  -  Af(x)) 

£  2 

A^f(x)  =  —  (A^f(x+6)  -  A^f(x)) 

L3 

where  the  scaling  factors  Ej,  &2>  e3  may  be  different  for  each  difference. 
This  method  would  permit  the  iteration  to  be  more  closely  tailored  to  the 
specific  polynomial  being  generated,  but  would  require  a  more  thorough 
and  time-consuming  analysis  in  order  to  select  the  proper  values. 
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A  further  variation  would  be  to  make  the  difference  6  different 

for  each  A^f,  computing  the  more  rapidly  changing  quantities  f,  Af 

more  frequently.  This  would  be  useful  in  the  hardware  implementation 

discussed  next,  as  it  would  allow  a  relatively  smooth  movement  of  the 

beam  between  a  set  of  coarsely  spaced  points  t.  on  the  curve,  the 

2  3  1 

coarse  points  determined  by  A,  A  ,  and  A  ,  with  f  itself  making  a 
smoothly  interpolated  transition. 

Roberts  and  Blatt  [4,  16]  have  designed  and  built  at  the  M.  I.  T. 
Lincoln  Laboratory  a  curve  generator  using  rational  parametric  quad¬ 
ratics.  It  uses  a  digital  counter  to  generate  a  variable  t,  0  <  t  <  1, 
and  then  uses  multiplying  digital  to  analog  converters  to  compute  the 
polynomials,  with  a  feedback  loop  to  perform  the  homogeneous  division. 
Figure  2.2.2  presents  a  suggestion  for  a  curve  generator  using  the 
iterative  approach  analyzed  here.  Conceptually,  it  would  be  somewhat 
better  than  Roberts'  scheme,  using  only  three  (optionally  four)  multi¬ 
plying  D /A  converters  with  only  one  converter  in  the  feedback  loop, 

whereas  Roberts  must  include  all  the  converters  a,  b,  c  making  up  a 

2  2 
sum  at  +bt+c,  as  well  as  the  converter  generating  t  and  t. 

The  additions  would  be  done  at  fixed  clock  intervals,  either 
serially  or  in  parallel,  depending  on  your  budget  and  desire  for  fast 
adders.  The  sign  of  the  feedback  amplifier  would  have  to  be  chosen  to 
make  the  loop  stable,  that  is,  opposite  to  the  sign  of  hz.  The  third 
column  is  optional  and  gives  an  intensity  signal  for  depth.  Digital 
scaling  by  e  would  be  performed  by  having  the  lower  inputs  to  each 
adder  shifted,  either  by  a  fixed  amount  or  selectable  from  a  set  of 
fixed  amounts,  say,  two  to  six  bits,  meaning  e  =  tt|-  .  Analog 
scaling  would  be  done  so  that  analog  overflow  would  indicate  either  a 
point  offscale  on  x  or  y  or  too  close  to  the  observer  in  z;  such  off- 
scale  portions  of  the  curve  would  be  blanked  out. 
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INTENSITY 


FIG-  2.2.2.  CURVE  GENERATOR  USING  MULTIPLYING 
D/A  CONVERTERS. 
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2.3  BASIC  CURVE  FORM 


N 

For  our  purposes  a  curve  in  R  is  a  vector-valued  function  from 

N  * 

the  real  line  R  to  the  real  vector  space  R  of  dimension  N.  That  is, 
a  curve  is  an  ordered  N-tuple  of  real  functions  (f^fg,  .  ..,fj^).  Further¬ 
more,  let  us  restrict  the  functions  f.  to  be  linear  combinations  of  some 
linearly  independent  set  of  functions,  <£  =  [  <i>^,  ^ ,  ...  ,  <j> ]  and  let  the 

functions  be  defined  over  the  domain  fO,  ll.  Let  f  (u)  =  S  A..6.(u) 

J  i=0  XJ  1 

where  the  {A^}  are  constants.  Then  a  curve  segment  can  be  represented 
by  a  unique  (M+l)  X  N  real  matrix  A.  A  point  V(u)  on  the  curve  associ¬ 
ated  with  a  parameter  u  is  given  by 

V(u)  =  [  $0(u),  <f> ^ (u) ,  ...  ,  <£m(u)]  A  =  <j>( u)A  = 

=  [f1(u),  f2(u),  ...  ,  fN(u)] ,  0  u  <  1. 


In  this  context  we  will  say  either  that  the  matrix  A  represents  the  curve, 
or,  more  pointedly,  that  the  matrix  A  ^s  the  curve. 

-►  N  L 

If  we  have  some  mapping  H:R  -*•  R  from  N  dimensional  into  L 

dimensional  space,  then  the  image  v(u)  =  H(V(u))  of  the  curve  V  under 

the  mapping  will  be  a  curve  v(u)  in  R  .  The  matrix  A  does  not  uniquely 

represent  the  curve  v(u).  We  can  define  the  equivalence  class  H[A]  of 

a  matrix  A  under  the  mapping  as  the  set  of  all  matrices  AM  such  that 

3("$(u)  A„)  =  3(<^(u)  A);  it  is  clear  that  all  the  curves  in  R  speci¬ 
al  T 

fied  by  the  A^eH[A]  will  be  mapped  into  the  same  curve  in  R  J  and 

there  will  thus  be  many  representations  of  that  curve.  Notice  that 

although  the  original  function  space  --  [f . ,  f2,  ...  ,  f^  ]  --  is  a  vector 

space  of  rank  not  exceeding  (M+l)  its  image  under  the  mapping  3  is  not 

necessarily  of  finite  rank  and  we  cannot  in  general  represent  the 

resulting  curve  in  R^  by  an  expression  of  the  form  W(u)  =  ^(u)B  where 

W(u)  e  R  ,  tf(u)  e  R  and  B  a  constant  matrix  of  dimension  P  X  L,  ^(u) 

a  fixed  set  of  basis  functions,  dependent  only  on  H  and  <f>.  The  mapping  3 


In  other  words,  we  only  consider  curves  in  a  parametric  (explicit), 
rather  than  implicit,  form. 
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will  of  course  map  the  curve  [f^,  fg,  ...  ,  f^  ]  into  a  curve  [g^,  gg,  •••  ,  ] 

where  g^  =  . ..,  f,^T],  where  the  {IT}  are  the  components  of  the 

vector  function  H. 

We  will  usually  choose  the  set  of  basis  functions  [<£n,  <£,,.]  to  be 

M  M-l  2  u  ivi 

the  polynomials  [u1  ,  u  ,  . ..,  u  ,  u,  l]  or  some  other  set  of  linearly 

independent  polynomials  of  degree  M.  The  functions  T  will  thus  be  the 

real  polynomials  of  degree  M  or  less.  We  will  also  choose  N  to  be 

dimension  3  or  4  and  L  to  be  dimension  2  or  3  depending  on  whether  we 

are  talking  about  two  or  three  dimensional  curves.  The  mapping  H 

will  be  that  used  in  making  the  transformation  from  homogeneous 

coordinates  to  ordinary  coordinates,  namely,  the  projection  obtained  by 

dividing  the  first  L  components  of  the  vector  V  (of  dimension  L+  1  =  N) 
s  t 

by  the  L+l  component  to  arrive  at  the  L  components  of  the  vector 
H(V).  We  thus  identify  the  curve  specified  by  the  form  <f>A  and  the 
resulting  vector  V(u)  as  homogeneous  coordinates  of  a  curve  v(u)  in 
dimension  one  less.  The  many-to-one  property  of  the  mapping  0[  is 
illustrated  by  the  fact  that  two  points  in  homogeneous  coordinates,  ? 
and  <5,  are  equivalent  if  and  only  if  ?  =  aC§  for  some  non-zero  a. 

Thus  we  talk  about  curves  of  the  form 


x 


M  .  M-l 

aMU  aM-lU 


,  M  .  ,  M-l 
dMU  +  dM-lU 


+  .  .  .  +  a^u  +  a 
+  .  .  .  +  d^u  +  d 


0 


0 


y 


,  M  .  ,  M-l 

bMU  +  bM-lU 


,  M  ,  ,  M-l 
dMU  +  dM-lU 


+  .  .  .  +  b^u  +  b^ 
+  .  .  .  +  d ^u  +  d^ 


and  optionally. 


z 


M  ,  M-l 

CMU  +  CM-1U 


,  M  .  ,  M-l 
dMU  +  dM-lU 


+  .  .  .  +  c  j  u  +  c 
+  .  .  .  +  d^u  +  d 


0 


0 


The  value  for  each  component  of  a  point  on  the  curve  is  the  ratio  of  two 

polynomials  in  the  parameter  u.  In  theory  and  in  practice  we  will  first 

N 

do  some  operations  on  the  curve  in  R  (or  upon  its  representation  as  a 
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matrix)  and  only  at  the  last  moment  before  display  perform  the  mapping 
N  L 

H  :  R  -*•  R  by  doing  the  division  to  remove  the  homogeneous  coordinate. 
In  most  cases  the  many-to-one  property  of  H  will  not  come  into  play 
explicitly  but  only  when  we  remember  that  the  points  V(u)  are  expressed 
in  homogeneous  coordinates. 

If  we  have  two  sets  of  basis  functions  for  the  original  function 
s pace  [ f j ,  •  •  •  ,  f-^ say  [$q>  0  j  >  •••  >  1  3-nd  [  ‘ 

always  write  the  curve  in  two  forms. 

V(u)  =  [0Q,  <t>v  ...,  <(>M] A  =  [¥Q,  tfj,  ...,  ^m]B,  B  =  TA  where  T  is 

the  change  of  basis  transformation,  [  \&q,  \k j,  . . .,  T  =  [  <£q ,  (j>^,  . ..,  • 

This  relation  will  prove  extremely  useful  in  the  practical  matter  of 
determining  the  matrix  A  to  represent  a  curve  under  the  basis 

=  [  </>0J  where  will  be  chosen  so  as  to  make  the  display 

easy  where  a  different  basis  might  be  more  convenient  for 

specifying  the  curve  or  performing  computations  upon  it. 


2.4  REPARAMETERIZATION  AND  SHAPE  INVARIANCE 


Supposing  we  have  a  curve  represented  by  the  matrix  A,  it  is 

natural  to  ask  what  matrices  B  represent  the  same  curve.  The  major 

reason  for  investigating  this  shape  invariance  is  to  enable  us  to  draw 

a  particular  curve  in  as  smooth  a  manner  as  possible,  while  using 

equally  spaced  {mj  .  The  initial  choice  of  a  matrix  to  represent  a 

curve  may  be  poorly  made,  as  when  the  velocity  vector  is  great  at 

C'U 

the  same  time  the  curvature  of  the  curve  is  great.  We  prefer  to  have 
the  velocity  slow  at  such  points  in  order  that  successive  points  in  time 
along  the  curve  well  represent  its  shape.  (See  Figure  2.4.1.) 


We  define  two  curves  in  R^  to  be  the  same  if  and  only  if  there 
exists  a  one-to-one  correspondence  between  the  points  on  the  curves 
such  that  for  each  point  "p(s)  on  curve  A  there  exists  a  point  q(u)  on 


We  normally  think  of  the  parameter  u  as  being  time,  and  derivatives 
with  respect  to  it  as  velocities,  for  reasons  to  become  obvious  later. 
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Fig.  2.4.1  Good  and  Bad  Parameterizations 
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curve  B  such  that  p(s)  =  q(u).  This  correspondence  between  the  points 
on  the  curves  induces  a  correspondence  between  the  parameters  on  the 
curves,  which  we  will  denote  by  a  function  s  such  that  s  =  s(u),  or 
p(s)  =  p(s(u))  =  q(u).  In  this  context  we  may  have  to  extend  the  domain 
of  one  or  the  other  of  the  curves  to  cover  the  entire  real  line, 

<  u  <  +°°,  <  s  <  +o°,  although  we  will  still  only  display  a  finite 

portion,  V(u)  =  <^(u)B,  and  0  5  u  2=  1.  (See  Figure  2.4.2.) 


In  the  notation  of  the  previous  section,  we  have 

rtr  r  M  M-l  1 1  a!  ft  r  r  M  M-l  .  -i  _-i 

H{  [  s  s  ...  s  1]A)  =H{[u  u  ...  u  1JB} 

where  H  is  the  operation  of  performing  the  homogeneous  division.  This 

equality  between  the  coordinates  on  the  curves  will  hold  at  the  point 

p(s(u))  =  q(u)  if  and  only  if  there  exists  a  non-zero  scalar  factor  a 

associated  with  the  point  p(s(u))  =  q(u)  and  hence  a  function  n(u)  of  the 

parameter  u  such  that 

/  \  r  ^1  M- 1  i  i  a  r  ^1  M- 1  i  i  -d 
q'(u)[s  s  ...  s  1 J A  =  [ u  u  ...ulJB. 

Write  the  above  equation  for  M+l  different  values  of  u,  and  get: 


cKu^) 

r  / 

*uo 

,  .M-l 

S(U0) 

.  .  s(u0) 

1 

^(Uj) 

,  ,M 

stUj) 

,  .M-l 

siuj) 

stu^ 

1 

^uM)_ 

,  •  ,M 

L  s(uM) 

/  \M-1 

S  UM 

s(uM) 

1_ 

u 


u. 


M 

0 

M 


u 


.M 

M 


u 


M-l 

0 

M-l 


u, 


u 


M-l 

M 


The  matrix  which  precedes  the  B  matrix  is  a  Vandermonde  matrix  [2], 
which  is  non-singular  since  the  {m}  are  all  different.  Note  that  the 
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s=s(  u ) 


P(s)  =  <£(s)  x  A 


~>  — > 

Q(u)  =  <^>  ( u )  xB  =  <£  ( u )  xS  xA 


Fig.  2.4.2  Shape  Invariant  Transformation 
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matrix  preceding  the  A  is  also  a  non-singular  Vandermonde  matrix,  as 
different  {m}  imply  different  {s^}. 

Solve  for  B: 


'  M 

uo 

1 

-1 

u(uQ) 

r  , 

s(u0) 

1 

M 

.  UM  '  • 

1 

_ 

Q'(um)_ 

/  nM 

LS(V  ■  ' 

1 

which  we  write  as  B  =  SA  for  some  as  yet  undetermined  square 
(M+l)  X  (M+l)  matrix  S.  Since  S  is  the  product  of  three  non-singular 
matrices,  it  is  non- singular . 


The  original  equation  becomes 


,  v  r  M  M-l 
nlu^s  s 


1  I  a  r  M  M  - 1 
1 J  A  =  [  u  u 


1  ]  SA 


which  holds  if: 


.  ,r  M  M-l 
a  (u)  [  s  s 


1  1  r  M[  M!  —  1  .  1  n 

1J  =  [u  u  ...  1JS 


If  the  rank  of  the  null  space  of  A  is  zero,  then  this  will  be  a  necessary 
condition  as  well.  We  will  not  pursue  a  further  characterization  of 
curves  in  terms  of  the  null  space  of  their  associated  matrix. 


Under  these  conditions  equating  appropriate  components  of  the 
vectors  in  the  above  equation  gives: 

(1)  for  the  last  (M+1)S^  component, 

0 


M 


a(u)  1  =  [u  ...  l]S 


0 
1  J 


=  P(u) 


P  is  a  polynomial  of  degree  not  greater  than  M.  If  P  =  u(u)  is  a 
constant,  then  s  =  u  and  S  =  kl  for  some  constant  k,  a  trivial  and 
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uninteresting  equality  of  curves.  Hence  we  assume  P  is  not  a  constant. 
(2)  for  the  component: 


M 


Q'(u)  S  =  [u  ...  l]S 


0 

1 

0 


=  Q(u) 


Q  is  a  polynomial  of  degree  not  exceeding  M  and  s  =  = 


(3)  for  the  first  component: 


/  \  r  iic 

cn (u)  s  =  [u  ...  1 JS 


1 

0 


=  R(u) 


R  is  a  polynomial  of  degree  not  exceeding  M  and  s  =  —  ^  =  p 


Q 

P 


M 


R 


The  expression  p  = 


Q 

P 


M 


has  at  most  M  poles  and  at  most  M 
zeros  (including  multiplicities)  since  the  polynomials  R  and  P  are  of 
degree  not  greater  than  M.  Thus  p  has  at  most  one  zero  and  one  pole 
and  we  have 

_  Q  _  au  +  b 
S  P  cu  +  d 

Substituting  in  the  original  equation  gives: 

M  /_  .  ,  .  M-l 

...  1 


,  .  /  au  +  b  \  /  au  +  b  \ 

“(u)Lvur+d)  IcuTd ) 


M 


=  [u1Vi  ...  1]  S  = 


(cu+d) 


M 


[(au+b)^1],  [(au+b)M  ^(cu+d)] ,  ...  ,  [(cu+d)^1] 


M 


Each  component  of  [ux  ...  l]S  is  a  polynomial  of  degree  not 
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exceeding  M  in  u.  The  components  on  the  right-hand  side  above  must 
thus  be  the  same  polynomials.  Hence  Q'(u)  =  e(cu+d)^,  where  e  is  a 
constant,  and  by  inspection  the  component  of  the  matrix  S  is  the 
coefficient  of  u^  *  in  the  expansion  of  e(au+b)^  ^(cu+d)^.  Bypassing 
the  algebra. 


J  k=0  J  /  / 


where  we  use  the  convention  that  (  =  0  when  j,  k  are  outside  of  the 

range  0  <  k  <  j. 


The  matrix  S  for 

M  =  2  is 

r  2 

a 

ac 

2  1 
c 

2ab 

be  +  ad 

2cd 

_b2 

bd 

d2_ 

and  the  matrix  for  M  =  3  is 

3 


3a2b 


3ab' 

,3 


2 

a  c 


2abc  +  a  d 


b  c+2abd 


b2d 


ac 

be2  +  2acd 
2bcd  +  ad2 
bd2 


3c2d 


3cd 

,3 


Typically,  we  will  want  to  use  the  reparameterization  matrix  S 
to  change  the  amount  (extent)  of  a  given  curve  that  is  drawn  and  to 
alter  the  rate  at  which  the  moving  vector  [u^  ...  l]SA  traverses  the 
curve.  These  variables  can  be  specified  by  setting 


1)  s(0)  =  a 

2)  s(l)  =  0 


3) 


s'(0) 

sTl) 


2 

r 


a  is  the  value  of  the  parameter  s  on  the  original  curve  A  at  which 
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the  new  curve  B  begins.  /3  is  the  value  of  the  parameter  s  on  A  at 

2 

which  the  curve  B  ends,  r  is  a  ratio  of  derivatives;  the  significance 
of  r  will  be  shown  shortly.  It  is  easily  verified  that  the  ratio 


s'tuj) 

s'(u2) 


cu2  +  d 

cUj  +  d 


5*  0 


a,  b. 


Substituting  the  above  conditions  into  s  = 
c,  d: 


we  solve  for 


a  =  j3  -  a  r 

b  =  ar  (B  -  £?r)u  +  ar 

or  S  =  -Jj-Z - r — T - 

.  ( 1  -  r)u  +  r 

c  =  1  -  r 

d  =  r 


using  the  arbitrary  condition  c  +  d=  1  to  remove  the  extra  degree  of 
freedom  introduced  by  the  homogeneous  form  of  the  function  s.  The 
parameter  e  will  allow  an  absolute  adjustment  of  the  homogeneous 
coordinate  system,  if  that  is  desired. 

The  velocity  vector  ^  at  the  beginning  (u  =  0,  s  -  a)  will  be  multi- 

plied  by  the  factor  —  and  at  the  end  (u  =  1,  s  =  jS)  by  (/3-a)r.  In  fact, 
in  general,  we  have 

dp~  _  ds  d£  _  (^-apr  dp~ 

du  "  du  ds  ■  [U.r)u+r]2  ds 

The  product  of  the  magnitudes  of  the  velocity  vectors  at  the  ends  of  the 

2 

new  curve  B  is  thus  (fi-a)  times  the  product  of  the  magnitudes  of  the 
velocity  vectors  at  the  corresponding  points  s  =  a,  s  =  j3  of  the  original 
curve  A.  This  product  thus  remains  constant  for  all  similar  curves 
covering  the  same  interval  ( a ,  /3)  of  A.  It  is  in  any  case  independent  of 
the  rate  parameter  r.  Notice  that  if  r  <  0,  the  new  curve  will  not  com¬ 
pletely  overlap  the  old  curve  in  the  range  0  ■*£  u  <  1,  0<s<l  although 
over  the  full  line  -°°  <  u  <  +°o,  -°o  <  s  <  +°o,  the  two  curves  will  be 
identical.  In  particular,  at  the  value  u  =  —  _  j ,  s(u)  =  +°°,  definitely 
outside  the  domain  of  the  original  curve  v(s). 

Thus,  for  example,  if  the  matrix  A  draws  half  a  circle,  the 
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matrix  SA  will  draw  the  other  half  in  the  same  direction  if  or  =  1,  /3  =  0, 
r  =  -1.  (See  Figure  2.4.3  for  an  example  of  the  use  of  this  continuation 
matrix.) 


For  M  =  2  this  matrix  S(l,  0,  -1)  is 

"1  2  4" 

-2  -3  -4 

1  1  1  _ 

and  for  M  =  3  it  is 

"1  2  4  8" 

-3  -3  -8  -12 

3  4  5  6 

_-l  -1  -1  -1_ 

The  matrix  S(l,  0,  1)  would  draw  the  same  curve  as  A  but  in  the 
opposite  direction.  This  points  out  the  important  fact  that  the  curves 
are  directed,  which  must  be  taken  into  consideration  when  they  are 
composed  into  surfaces,  where  the  directions  of  opposite  boundary 
curves  makes  a  big  difference  on  the  surface  between  them. 


If  we  wish  only  to  change  the  parameterization  of  a  curve  without 
changing  its  extent,  we  have  a  =  0,  /3  =  1,  and  the  reparameterization 
matrix  for  M  =  2  is 


S(r)  = 


1 

0 

0 


1  - r  (l-r)^ 
r  2r(l-r) 


and  for  M  =  3 


1 

(1-r) 

(1-r)2 

(1-r)3 

S(r)  = 

0  r 

2(l-r)r 

9 

3(  l-r)2r 

o 

0 

0 

r 

3(l-r)r 

L0 

0 

0 

3 

r 
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Fig.  2.4.3  Continuation  of  a  Curve 
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This  matrix  can  be  used  to  adjust  the  homogeneous  coordinate 
representation  of  a  curve  in  an  essentially  arbitrary  manner,  as  we 
now  demonstrate.  Let  us  assume  a  curve  of  the  form 


P(u)  =  [u^u^ul] 


a 

b 

c 

d 


r  3  2  i  1  A 

=  [u  u  u  1 J A 


that  is,  P^(0)  =  d,  P^(l)  =  a  +  b  +  c  +  d.  Applying  our  matrix  to  the 
above,  we  obtain 

Q(u)  =  [u3  u3  u  l]SA  =  [u^u^ul] 

where 

a  =  a  +  b(l-r)  +  c(l-r)3  +  d(l-r)3 
0  = 

7  = 

6  = 

and  thus 

Qh(0)  =  6  =  dr°  =  r°Ph(0) 

Qh(l)  =  a  +  f 3  +  7+  6  =  a  +  b  +  c  +  d  =  Ph(l). 

If  we  also  apply  a  rescaling  of  the  homogeneous  coordinates  by 
R(u)  =  sQ(u),  we  have 

Rh(0)  =  sr3Ph(0) 

Rh(l)  =  sPh(l) 

and  we  can  make  R^(0)  and  R,  (1)  completely  arbitrary  by  setting 


b(r)  +  2c(l-r)r  +  3d(l-r)3r 


cr 


+  3d(l-r)r‘ 

,  3 
dr 


...  a 
•  .  .  0 
•  •  •  7 
.  .  .  6 
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h( 0)  =  Rh(0)  -  sr3Ph(0) 
h(l)  =  Rh(l)  =  sPh(l) 


and  solving  for 


-h(0)Ph(l)- 
_h(l)Ph(0)  . 


1/3 


and 


s 


h(  1) 

ph(D 


so  long  as  none  of  the  quantities  are  zero.  In  summary,  given  an  arbi- 

3  2 

trary  curve  P  =  [u  u  u  l]  A,  we  obtain  the  same  curve  in  the  form 
3  2 

R  =  [u  u  u  l]sS(r)A  where  we  have  completely  specified  the  homo¬ 
geneous  coordinates  as  R^(0)  =  h(0),  R^(l)  =  h(l).  In  other  words,  there 
are  a  doubly  infinite  number  of  matrix  representations  for  a  given  curve 
segment. 

If  we  have  a  curve  v(u)  =  15  [^A^,]  represented  by  a  matrix  A,p  for 
some  other  basis  function  ip  =  [u  u  ...  l]T  then  the  reparameteri¬ 
zation  matrix  S™  for  this  basis  will  be  given  by  the  similarity  transfor- 

- 1  1 

mation  S,p=T  ST  and  the  curve  matrix  after  reparameterization  will  be 

t_1stat. 


2.5  ENDPOINT  DERIVATIVE  FORM 

In  the  previous  subsections  we  have  been  talking  about  curves  of 
arbitrary  degree;  in  this  subsection  we  will  be  concerned  only  with 
curves  of  degree  m=3.  Let  hv  =  V  denote  the  homogeneous  coordinates 
of  a  point  in  R^;  that  is,  let  h  =  V^,  v  =  where  v  is  a  vector 

of  the  form  [viv2*  •  •  vn_  i  l]  Then  the  -^equation  of  a  curve  is  given 
parametrically  by 

V(u)  =  h(u)v(u)  =  <f>(  u)A 

where  <j>( u)  are  appropriate  M-dimensional  basis  functions.  The  curve 
in  R  1  =R  J  is  obtained  by  taking  the  first  L  components  of  the  vector 
v(u);  that  is,  by  dividing  out  the  homogeneous  coordinate  h(u)  from  the 
vector  h(u)v(u)  and  dropping  the  last  component. 
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Let 


V(0)  =  V 
V(l)  =  V, 


0 


vmo)  =  dyu> 

du 

v'(D  =  v; 


=  V 


u=0 


then  we  have 


"  V 

"?(0)“ 

V1 

?(1) 

% 

?'(o) 

_n. 

£'(1) 

One  can  show  that  the  square  matrix 

?(0>  ' 

?(1) 

?'(0) 

?'(1)  J 

is  non-singular  since  the  basis  functions  $”(u) 
linearly  independent  polynomials  of  degree  3. 


M  = 


and  then 


A  =  M 


0(0) 

0(1) 

?'(0) 

?'(D- 

vo 

V, 


vf 

0 


-1 


are  by  definition  four 
Define 
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We  now  write 


vo  =  hovo 


v!  =  Vi 

T^,  _  d(hv ) 
0  du 


u=0 


=  hf  v  +  h  v  ' 
0  0  0  0 


V'  =  h' v  +  h  v  ' 
1  11  11 


where  v^,  are  the  derivatives  with  respect  to  the  parameter  u  and 
where  h^,  h^,  h|()  h^  are  meaningful  only  in  homogeneous  coordinates. 
Notice  that  the  last  component  of  v^  and  v^  is  now  0.  Thus 


A  =  M 


h0 


0 

hl 

0 

h' 


0 

0 

h. 


1 


1 


The  curve  is  given  by 


1(u)A  =  1(u)M 


h. 


h' 


h. 


hi 


h. 


h. 


1 


We  can  perform  the  multiplication  <jT(u)M  separately  to  define  a  new  set 
of  basis  polynomials  [F^F  G^Gj]  =  <f>  M  with  transformation  of  basis 
matrix  M  itself.  We  observe  that 


_F0(0) 

F^O) 

G0(0) 

G^of 

’1(0) 

1 

0 

0 

0 

F0(l) 

f^d 

Go(1) 

Gjd) 

?d) 

II 

§ 

i-H 

1 

§ 

ii 

§ 

0 

1 

0 

0 

FJ}(0) 

Fj(0) 

Gyo) 

Gj(0) 

?'<0) 

0 

0 

1 

0 

( - 

1—* 

Fjd) 

G’(l) 

G'jd) 

1 

1 

Pd) 

L. 

0 

0 

0 

1 

This  notation  is  from  Coons  [7]. 
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With  these  basis  functions,  a  curve  can  always  be  represented  by 


a  matrix  of  the  form 


AM 


h. 


hr 


h. 


h. 


h. 


Ll 


since  it  is  given  by  the  equation 


V(u)  =  [FqFjGqGJ 


h. 


h. 


h. 


h. 


h. 


h. 


We  will  sometimes  write  this  endpoint-derivative  form  as 


AA(u)  =  [FqF^GJ 


(hv )( 
(hv). 
(hv)j 
(hv) 


1 J 


It  will  also  prove  useful  to  rewrite  it  in  the  form 


V(u)  =  [h0h  HJliJ] 


i  j 


or  in  standard  form  as 
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V(u)  =  [(FQh0  +  GQh^),  (F1h1  +  Gjh^),  (GQh0),  (G^J] 


V1 
v~0  ! 


where  now  the  matrix 


1  i 


A  = 


is  especially  simple,  although  for  a  very  special  basis  function  con¬ 
structed  specifically  for  the  particular  set  of  homogeneous  coordinates 
h^,  h^,  hp,  h!^.  This  form  shows  that  a  point  on  the  curve  V(u)  is 
always  a  linear  combination  of  the  vectors  Vq,  v^,  Vq,  Vy 

For  these  basis  functions  the  reparameterization  matrix  becomes 


S(r) 


r3  0  0  0 

0  10  0 
3(l-r)r2  0  r2  0 

0  3(l-r)  0  r 


using  the  similarity  transform  M  ^SM. 

Notice  that  this  result  gives  a  particularly  simple  answer  to  the 
question  of  determining  whether  two  matrices  represent  the  same 
curve,  namely,  if  in  endpoint  derivative  form  we  have 

3(0)  =  O'r3^(0) 

3(1)  =  aM  1) 

3'(0)  =  a  3(  l-r)r2  ^(0)  +  qt2^,(0) 

3'(1)  =  a  3(l-r)I(l)  +  or^'(l) 
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then  and  only  then  will  the  matrices  A  and  B  represent  the  same  curve. 
a  and  r  are  arbitrary  constants. 


2.6  CURVE  SPECIFICATION 


In  the  previous  subsection  we  have  characterized  a  rational  para¬ 
metric  cubic  polynomial  curve  in  terms  of  endpoints  --  v^,  v^  --  and 
tangent  vectors  at  the  endpoints  —  v^,  v^  --  with  four  remaining  degrees 
of  freedom,  h^,  h^,  h^,  h^.  In  this  subsection  we  investigate  several 
ways  of  computing  the  numbers  h^,  h^,  h^,  and  h^  in  order  to  completely 
specify  the  curve  from  geometric  consideration.  (See  Figure  2.6.1.) 
Although  some  or  all  of  these  results  may  be  inappropriate  for  a  par¬ 
ticular  application,  the  techniques  used  indicate  the  generality  of  the 
curve  formulation  and  what  we  believe  to  be  the  proper  way  to  attack 
the  problem. 


Let  us  require  the  curve  to  pass  through  the  point  v^  at  the  value 
u^  of  the  parameter.  We  have 


7C  =  lhohihohi) 


F0(uc> 


G0(uc) 


F1<V 


°l(uc) 


°0(uc> 

0 

0 

0 


Gl<uc> 

0 

0 


I  V. 


v . 


=  hFV 


v 


1  i 


We  can  solve  =  hFV  for  hF  (and  hence  for  h,  since  F  is  non¬ 
singular  for  any  0  <  u  <  1)  if  and  only  if  v  belongs  to  the  range  of  V. 

3  c  c 

In  R°  the  geometrical  meaning  of  this  condition  is  that  if  the  curve  is 

planar,  v  must  be  in  the  same  plane,  or,  if  the  curve  is  linear,  v 
must  be  on  the  same  line.  If  the  curve  is  really  three-dimensional, 
then  V  will  have  full  rank  and  we  can  solve  directly  for 

h  =  v  V'V1  =  v  (FV)"1 
c  c 


2  — 

In  R  the  equation  v  = 

equations  and  four  unknowns. 


hFV  is  indeterminate  since  we  have  three 
If  v£  is  in  the  range  of  V,  we  can  impose 
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/ 


/ 


Fig.  2.6.1  Specifying  the  Curve 
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an  additional  condition.  ’  One  convenient  such  condition  is  to  specify  that 
— ►  x* 

the  slope  at  v  be  given  by  the  ratio  _ c  _  dx  .  Taking  derivatives  with 

y'c  "  dy 

respect  to  u,  we  have 


[hohihohi^ 


Fo 


Go 


G! 


G0 


G! 


1 


d(h  v  ) 
c  c 

du 


u=u 


=  hf  v  +  h  v  '  =  h'  v  +  v ' 
c  c  c  c  c  c  c 

because  we  have  let  h^  =  1.  (v£  =  [x^y^l].)  Since  the  third  component 

of  vf,  vT  vJ  is  0  we  can  solve  for 
c’  0'  1 


h0  =  lhohihohi) 


Ff 

0 


g: 


g: 


and  thus 


"Fo 

Go 

vo 

'Fo~ 

[h0hlhbhi) 

Fi 

G0 

Gi 

vo 

-  (hlhlhohil 

Fi 

Go 

.  Gi 

,GiJ 

v  +  v' 
c  c 


or,  rearranging, 


vc  will  not  be  in  the  range  of  V  if  and  only  if  v^,  v^  and  vq  “  are 

collinear  and  v  -  v,  is  not  on  the  same  line.  This  means  that  the 
cl  _ 

tangents  at  each  endpoint  point  to  the  other  end  and  v£  is  not  on  the 
line  between  the  two  endpoints. 
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[ hohihohi  1 


F' 

0 


Denote  the  left  side  of  the  above  by  Q.  Since  we  are  only  con¬ 
cerned  with  the  slope  at  v  ,  we  can  extract  the  first  and  second 
components  by  a  post-multiplication  and  obtain  the  ratio 


Q 


"  l" 

o 

0 


Q 


"O' 

l 


o 


whence 


Q 


y 


-X' 


0 


0 


If  we  adjoin  this  single  equation  to  the  original  equation  for  the  curve 
passing  through  the  point  v  ,  we  have 


Ff 

0 


If  the  matrix  is  invertible,  we  can  solve  for  h^,  h. ,  h^,  h^  and 

—i ►  x  f 

determine  a  curve  passing  a  specified  point  v  with  a  given  slope  c 

c  • 

y 

J  c 

— ►  — — ►  — ► 

with  given  endpoints  v^,  v^  and  tangent  vectors  v'  v^  at  those 

endpoints  As  with  the  three-dimensional  case  above,  we  could  charac¬ 
terize  this  problem  in  somewhat  more  generality  in  terms  of  the  range 
of  the  appropriate  matrices,  but  we  will  not  do  so  since  the  conditions 
are  not  simple  enough  to  be  interesting. 

We  could  in  the  previous  examples  have  just  as  easily  required 


Fjo) 
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the  curve  to  pass  through  four  non-planar  points  v at  dif¬ 
ferent  values  a,  (3,  y,  6  of  the  parameter.  We  would  have 


T(o') 

h  v 

A  = 

O'  O' 

h/,s 

^(7) 

h  v 

7  7 

?(6) 

—  — 

Av 

^(u) 

-1 

!>* 

l _ 

Ve 

<£(7) 

h  v 

7  7 

?(6) 

Jv«_ 

Then  let  us  ask  that  the  curve  pass  through  vp  at  u  =  uc: 


^(uc)M 


h  v 
a  a 

h  v 
7  7 

h/«j 


v 

c 


where  now  we  have  let 


M  = 


^(a)  1 

<£(/3) 

<£(7) 

?(6) 


If  we  define  ^(uc)M  =  [F^F-^FgF^],  we  have  after  rearranging 
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which  can  be  solved  directly  for  [ h^h^h^h^ ]  if  the  matrix  is  non¬ 
singular,  that  is,  if  the  four  points  are  non-planar.  In  this  way  we 
have  asked  the  three-dimensional  curve  to  pass  through  five  specific 
points  at  five  specific  values  of  the  parameter.  Similar  results  could 
be  derived  in  two  dimensions. 


SECTION  III 


SURFACES 


3.1  INTRODUCTION 

Having  more  or  less  exhausted  the  interesting  properties  of 
rational  parametric  cubic  curves,  we  turn  to  a  consideration  of  a  class 
of  similarly  defined  three-dimensional  surfaces.  After  presenting  the 
basic  definition  of  this  class  of  rational  bi-cubic  surfaces,  this  section 
presents  several  different  but  related  characterizations  of  them.  In 
general,  we  define  a  surface  by  means  ofa4X4X4  =  64  element 
tensor.  The  difficult  problem  is  to  determine  the  numbers  to  place  in 
the  tensor  in  order  to  display  a  particular  desired  surface.  We  will 
talk  about  characterizing  the  surface  in  terms  of  its  boundary  curves, 
tangent  vectors,  normal  vectors,  and  points  through  which  the  surface 
is  to  pass.  We  do  not  speak  about  a  complete  surface,  but  rather  about 
patches  which  may  be  combined  to  form  a  more  complete  surface; 
questions  of  continuity  arising  in  this  connection  are  also  discussed. 

3.2  SURFACE  EQUATION 

Let  be  a  4  X  4  X  4  array,  then  the  components  of  the  homo¬ 
geneous  coordinates  of  a  point  on  a  surface  are  given  by  a  vector 
function  ?(u,  v)  of  two  parameters,  the  two  degrees  of  freedom  on  the 
surface,  as 

P.  =  u*B.  .  v^  ;  0  <  u  <  1,  0  <  v  <  1. 

k  ljk  '  * 

This  expression  may  be  written  in  matrix  notation  as 


Phx  =  lu3u2ul)Bhx 
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p 


hy 


[u3u2ul]Bhy 


3 

v 

2 

v 

v 


1  J 


Phz  =  [  et  cetera] 


P^  =  [  et  cetera] 


These  expressions  may  also  be  written  without  subscripts  in  the  follow¬ 
ing  fashion  to  indicate  all  four  equations,  the  implication  being  that  B 
is  a  matrix  whose  components  are  vectors: 


?(u,  v)  =  [u3u2  u  l]  B 


3 

v 

2 

v 

v 


1 


We  could  also  define  a  time-varying  surface  as 


Pk(ujVjt)  =  u1Bijk(t)v;j 

where  B. (t)  =  t  T...  this  surface  would  change  smoothly  from 
ijk 

Bijk^O)  =  ^ijkO  *°  ®ijk^  =  ^  ^ijkj?  aS  *  Ganges  from  0  to  1.  Even 
more  simply,  we  could  define 


Pk(u,  v,  t)  =  u1[tB.jk+(l-t)Bf;.k]vJ* 

in  which  the  surface  would  change  "linearly"  from  the  surface  B  to  the 
surface  B'.  Further  generalizations  are  obvious. 

Having  defined  the  class  of  surfaces  to  be  considered,  the 
remainder  of  this  section  is  concerned  with  determining  the  64  numbers 
in  the  tensor  that  defines  the  surface,  subject  to  various  constraints. 
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3.3  ENDPOINT  DERIVATIVE  FORM 


Borrowing  Coons'  notation  [7],  let  F^,  F^,  Gq,  be  functions 
of  one  variable  such  that 


F^j)  =  ;  G.(j)  =  0 
F'(j)  =  0  ;  G'(j)  =  fij 


Then  we  can  define  a  surface  as 


i,  j  =  0,  1 


?(u,  v)  =  [  Fq(u),  F^u),  Gq(u),  G^u)]  A 


F0(v) 

Fj(v) 

G0(v) 

Gj(v) 


Because  of  the  definition  of  the  F^,  G^  functions,  the  components 
of  A  relate  directly  to  various  partial  derivatives  of  the  surface 
function  ?(u,  v)  as  indicated  below: 


00 

01 

00 

V 

01 

V 

10 

11 

10 

V 

11 

V 

00 

u 

01 

u 

00 

uv 

01 

uv 

10 

u 

11 

u 

10 

uv 

11 

uv 

where  the  notation,  from  Coons,  is  exemplified  by 


1 

uv 


9?(u,v) 

3u3v 


u=l 

v=l 


01 

u 


3?(u,  v) 

3u 


u=0 


v=l 


Note  well  that  these  derivatives  refer  to  each  component  of  the  four¬ 
dimensional  homogeneous  coordinate  vector;  for  example, 

^  _  3(h(u,  v)p(u,  v))  _  3h  -»  .  3p~ 

u  3u  3u  p  3u 

and  we  must  solve  for 
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9 £  = 

3u 


h 


where  p  =  [xyz  1], 


12  = 

3u 


[  x'  y'  z'  0] 


to  obtain  the  partial  derivatives  of  the  three-dimensional  components, 
x(u,  v),  y(u,  v)  and  z(u,  v). 


If  we  insist  on  F^,  being  cubic  polynomials,  the  conditions 
define  them  completely  as 


3  2 


[  F^ujF^uJG^uJG^u)]  =  [u  u“ul] 


2  -2 
-3  3 

0  0 
1  0 


1  1 

-2  -1 

1  0 

0  0 


The  4X4  matrix  above  occurs  often  and  has  become  identified  with  the 
symbol  M,  for  "magic"  matrix.  The  tensor  A  in  the  definition  of  the 
surface  is  thus  the  representation  of  a  surface  using  the  basis  functions 
Fq,  Fj,  Gq,  and  Gr 


We  can  convert  between  the  two  forms  of  representation  thus  far 
introduced  --  the  polynomial  basis  and  the  endpoint  derivative  basis  -- 
by  means  of  the  identities 


B  =  MAMT 
A  =  M-1BM 
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since  [FqFjGqG  ]  =  [u  u  ul]M.  In  these  cases  the  matrix  multipli¬ 
cation  implied  in  the  equation  is  that  obtained  by  taking  each  slice  of 

T 

the  tensors  separately,  as  B,  =  MA^^M  ,  or,  written  in  indicial  form. 


B...  =  M.„  A.  ,  M  . 
ljk  ijk  imk  mj 

A..,  =  (M-1).  „  B„  .  (M_1)  . 

i]k  l  St  fmk  mj 


3.4  CURVES  ASSOCIATED  WITH  THE  SURFACE 

By  fixing  one  parameter  a  curve  in  the  surface  is  traced  by  letting 
the  other  parameter  range  over  the  interval  [  0,  l]  .  For  example,  let 
v  =  b,  then 
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P(u,  b)  =  [uVul]B 


b 

b 


=  [u3u2ul]A. 


ub 


where 


A  ,  =  B 
ub 


b‘ 

b 

1 


is  the  matrix  representing  the  curve  P 
cubic. 


v=b 


as  a  rational  parametric 


The  curves  P(0,  v),  f(l,v),  f(u,  0),  P(u,  1)  are  the  boundary  curves 
of  the  surface;  in  the  shorthand  notation  they  are  denoted  by  Ov,  lv,  uO,  ul. 

We  will  denote  the  matrix  describing  a  curve  for  u  =  constant  by 

a  /r  3  2  a,  ™T 
A  =  ([  a  a  a  1]  B) 
av  1 

so  that  the  curve  is  given  by  an  expression  of  the  standard  form, 

P(a,  v)  =  [vVvl)A 

For  computer  display  of  the  surface  we  generate  the  family  of 
curves 


P  ,  ,  b  =  —  ,  i  =  0,  1,  .  .  .  ,  m 

ub'  m 


or  the  family 


P  ,  a  =  i-  .  i  =  0,  1,  .  .  .  ,  n 
5 

ub 


av'  n 

or  both.  The  continuous  curve  PiiK  will  be  approximated  by  a  series  of 
chords 


P  ,P  ,  ,  u.  =  -r  ,  i  =  0,  1,  ...  ,  k-1 

u.b  u. ,  .b’  i  k 
i  l+l 
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and  P  by  the  chords 
av  J 


av .  av . ,  , 
J  J+l 


V* 


j  =  0,  1, 


i- 1 


3.5  ITERATIVE  DISPLAY  OF  SURFACE 


The  successive  chords  of  a  curve  and  the  successive  curves  of  a 
family  on  the  surface  are  computed  using  one  of  the  finite  difference 
schemes  of  subsection  2.2.  If  we  use  forward  differences,  we  thus 
compute 


Pk<m«,nT>  -(“)  e1CijkKJ  (") 


for  an  appropriate  tensor  C  dependent  on  the  surface  and  the  particular 
family  of  curves  to  be  generated.  The  multiplication  by  e*  is 

equivalent  to  multiplication  by  the  matrix 


1  e 

0  1 

0  0 

_  0  0 


0 

e 


0 

0 


m 


1  e 

0  1. 


which  is  achieved  by  simple  shifts  and  additions,  and  similarly  for 

.  The  summations  over  i  and  j  are  commutative,  so  if  we  are 
drawing  the  curves  ?(u,  ny)  we  first  compute  a  4  X  4  matrix 
cjj^  =  C.  for  each  curve  and  then  iterate  upon  it  to  compute  the 

successive  points  P^(m6,  n7)  =  ^  ^  J  e  on  the  curve  We  actually 


compute  through  the  difference  equation 


cn+1  =  cn  +  KCn 

ijk  Sjk  +  KSj+lk 


where  and  =  c!1^-  Similarly,  to  compute  P^(m6,  ny), 

we  compute 

^nm+l  _  „nm  .  „nm 
L,  —  L . ,  +  tL 

lk  lk  l+lk 

where  Cn.^  =  Cn  and  P.  (m6,n7)  =  C^P1. 

i]  ij  k  ’  "  Ok 
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The  initial  tensor  of  forward  differences  C...  is  obtained  from 

ijk 

the  tensor  representing  the  surface  by  the  relation 


=  Du<6,  e) 


^imk 


D. 

jm 


(7,  k) 


where  the  D  matrices  are  obtained  from  the  second  expression  on 
page  26. 

To  draw  the  other  set  of  curves,  P(m6,  v),  we  first  compute 
C“  =  e^ijk  an<^  derate  to  compute  P^(m6,  ny)  =  . 

The  algorithms  for  using  backwards  differences  are  similar  and 
readily  determined  from  the  above  and  from  subsection  2.2. 


3.6  TRANSFORMATIONS  OF  THE  SURFACE 

In  the  previous  subsections  the  coordinates  of  a  point  ?(u,  v)  on 
the  surface  have  been  defined  by  expressions  of  the  form 


Pk(u,v) 


U.(u)S. .,  V.(v)  where  U  and  V  are  vector-valued  functions  of 
l  ljk  y 


the  parameters.  Let  T  be  a  projective  transformation  to  be  applied 
to  the  surface  at  each  point  ?  and  let  P'  be  the  transformed  point. 
Then  we  have 

pi  -  PkTW  ■  <Ui<“>SijkVj<V>)Tki  " 

-  ui(u><sijkTki)vj(v)  - 

-  ui(u)S!jAj(v) 

where  the  surface  tensor  S^k  is  replaced  by  a  transformed  tensor 


-  ‘-'jjk--1-  kjg-  This  derived  tensor,  SL^,  will  generate  the  trans- 


SC.,  =  S. .,  T, 

lji  ijk  1 

formed  surface,  without  the  need  to  transform  each  point  to  be  dis¬ 
played.  In  other  words,  the  sixteen  vectors  ST  =  j  1  ^ij 2"^i j 3 ^ 

transform  in  the  same  fashion  as  points  on  the  surface  and  can  be 
regarded  as  a  set  of  sixteen  points  which  define  the  surface.  Notice, 
of  course,  that  any  transformation  --  including  the  identity  trans¬ 
formation  oT  --  of  the  surface  must  be  applied  identically  to  all 
sixteen  £>...  These  properties  are  obviously  what  allow  us  the  right 
to  call  the  array  a  tensor. 
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Notice  in  particular  that  the  tensor  may  be  the  set  of 
differences  needed  to  generate  the  surface  through  finite  difference 
iterations  Transformations  of  the  surface  can  thus  be  made  after 
the  display  representation  has  been  chosen  and  need  not  be  applied 


to  the  abstract  representation  B. 


ijk' 


3.7  REPARAMETERIZATION 

As  with  a  curve  matrix,  the  tensor  describing  a  surface  can  be 
reparameterized  to 

1.  change  the  rate  at  which  the  surface  is  traversed  by  the 
parameters 

and  to 

2.  display  other  (smaller  or  larger)  portions  of  the  surface  with 
the  same  parameter  range,  0  5  u  i  1,  0  i  v  S  1. 

In  the  case  of  surfaces,  1  is  particularly  important  for  such  a 
reparameterization  changes  the  appearance  of  the  curvilinear  net  used 
to  display  the  surface,  by  changing  the  density  of  contour  lines  non- 
uniformly,  whereas  with  curves  it  affects  only  the  accuracy  with  which 
the  chords  approximate  the  curve. 

3  2 

For  each  of  the  basis  functions  [u  u  ul],  [F  F  GnG-]  or 

u  u  u  2  u  3  v  i  u  j. 

[<0)(l)6(2)6  (  ^) 6  ],  there  exists  a  reparameterization  matrix  S(ar,  /3,  r), 

a  function  of  three  parameters  a,  j3,  and  r  which  will 

1 .  map  u  =  0  to  u  =  a 

2.  map  u  =  1  to  u  =  j3 

3.  change  the  parameter  rate  by  a  factor  r. 

These  matrices  are  obtained  as  described  in  subsection  2.4. 

Two  of  these  reparameterization  matrices,  S(ar,  |3,  r)  and  S(y,  6,  s) 
can  be  applied  to  a  surface  tensor  T  to  compute  a  new  surface  tensor 
T'  as 

Tijk  =  Si/a>  P’  r)  Tj?mkSjm(7j  6j  s) 
in  which  the  parameterization  in  both  u  and  v  has  been  altered. 
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A  shape -invariant  transformation  of  the  form  S(0,  l,r)  combined 
with  the  identity  a  I  has  the  property  of  allowing  us  to  adjust  the  fourth 
homogeneous  coordinate  h  at  three  of  the  corner  points  at  will;  in 
particular,  we  could  arbitrarily  assign  n  ^  =  hi  =  h^  =  1  without  any 
loss  of  generality  in  the  class  of  surfaces  represented.  Correspond¬ 
ingly,  if  we  are  given  a  surface  constrained  in  this  fashion  we  can 
reparameterize  it  so  as  to  arrive  at  a  more  satisfactory  rate  of  display. 


To  make  this  arbitrary  assignment  of  values  to  the  homogeneous 

coordinates,  recall  that  by  the  results  of  subsection  2.4  we  can  adjust 

the  homogeneous  coordinates  at  the  ends  of  a  curve  at  will.  Suppose  we 

have  a  tensor  A. .,  representing  the  surface  A.  Applying  a  reparameter- 
11K 

ization  matrix  S(r)  to  u  we  obtain  B.^k  =  S-^rjA^ ^  in  which  at  the 
corners  we  have 


B°°  .  r3A00 
h  h 

B?1  =  rV1 
h  h 

B^°  =  A?;0 
h  h 

b*1  =  a“ 

h  h 

Then  applying  a  reparameterization  S(s)  to  v,  we  obtain 


C...  =  B.  ,  S.  (s)  =  S.„(r)A„  .  S.  (s) 

ljk  lmk  jm  ii  Umk  jm 


in  which 


c“°  =  s3B°°  =  s3r3A°° 
h  h  h 

C31  =  B31  =  r3A31 
h  h  h 

c}°  =  S3b}°  =  s3a1° 

h  h  h 

C?;1  =  B* 1  =  A^1 
h  h  h 


and  lastly  applying  a  complete  rescaling  D 


ijk 


we  have 
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D 

D 

D 

D 


00 

h 

01 

h 

10 

h 

11 

h 


3  3.00 
=  0-3  r 

3  A  01 
=  ar  A, 
h 


3  A  1  0 


a  A 


11 

h 


If  none  of  the  or  A^  are  zero,  we  can  obviously  pick  a,  s,  r 
in  order  to  assign  arbitrary  values  to  the  coordinate  h  at  three  of  the 
corners;  the  value  of  h  at  the  fourth  corner  will  then  be  determined. 
There  are  thus  a  triply  infinite  number  of  tensor  representations  of  a 
surface  patch. 


3.8  CONSTRUCTING  THE  SURFACE 

The  preceding  subsections  describe  some  of  the  properties  of  a 
given  surface  and  are  thus  analytical;  the  usual  question,  however,  is 
the  constructive  problem  of  generating  a  surface  from  external  con¬ 
ditions.  This  and  the  following  two  subsections  sketch  some  partial 
solutions  that  are  being  investigated 

Suppose  we  are  given  the  four  boundary  curves  of  a  surface, 
specified  by  matrices  A^v,  A*v,  AuU,  A11'*'  in  endpoint  derivative 
form.  By  the  geometrical  closure  of  the  curve  segments  bounding  the 
surface,  these  matrices  are  guaranteed  to  represent  the  same  corner 
points  P(0,  0),  P(0,  1),  P(l,  0)  and  P(l,  1)  at  the  intersections  of  the 
generated  boundary  curves.  This  representation  is,  however,  in 
homogeneous  coordinates  and  we  have  no  guarantee  that  the  four¬ 
dimensional  vectors  given  by  the  matrices  are  the  same.  (In  other 
words,  ^V(v=0)  is  proportional  to  but  not  necessarily  equal  to 
^  U^(u=0).)  By  reparameterizing  each  of  the  boundary  curves  sepa¬ 
rately,  we  can  assign  absolute  values  to  the  homogeneous  coordinate  h 
at  each  of  the  corners,  adjusting  the  curve  matrices  A^v,  A^v,  Au^, 
Au^  so  that  the  endpoints  will  now  be  represented  by  identical  vectors. 

That  is,  we  adjust 
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Au0 

so  that 

Au0  /  n,  ,00 
Ah  (u=0)  =  h  , 

«u0,  1N  ,10 
Ah  (u=l)  =  h 

Aul 

so  that 

Aul /  m  ,01 
Ah  (u=0)  =  h  , 

. ul ,  - ,  ,11 

A,  (u=l)  =  h 
h 

A0v 

so  that 

»0v,  n\ 

Ah  ^v=0'  h  ' 

A°v(v=l)  =  h01 

Alv 

so  that 

a lv,  ,  ,10 

A^  (v=0)  =  h  , 

A^v(v=l)  =  h11 

In  assigning  the  four  values  of  h,  only  one  degree  of  freedom  will 
affect  the  shape  of  the  surface.  The  other  three  will  be  used  in  speci¬ 
fying  the  parameterizations  in  u  and  v  and  in  assigning  an  absolute 
homogeneous  scale  to  the  surface.  We  could  thus  arbitrarily  assign 
n  ^  =  n  ^  =  h^  =  1  and  solve  for  n  ^  according  to  some  additional 
criterion. 


At  this  point  we  have  filled  in  the  tensor  as  shown  in  Figure  3.8.1. 

(For  example,  S00Q  =  P°°  =  A“°;  S3q2  =  <P°\y  -  A®J.)  The  remain- 

ing  sixteen  components  can  be  determined  from  any  sixteen  independent 

conditions.  These  sixteen  as  yet  unspecified  components  have  come  to 

be  called  the  twist  partition  of  the  tensor  since  they  refer  to  mixed 

partial  derivatives,  •  A  possible  set  of  conditions  would  be  to 

specify  fT  at  each  of  the  four  corners  and  to  specify  a  point  p  through 

which  the  surface  must  pass  at  some  specified  value  of  the  parameters, 

say  u  =  -jj,  v  =  i .  With  these  conditions,  we  can  solve  for  (P^)uv  at  the 

four  corners  with  one  degree  of  freedom  left,  say  to  specify  the  h  value 

reached  at  the  center  points  ?  or  some  other  appropriate  geometrical 

s 

constraint. 


Instead  of  specifying  p^v  directly,  one  could  specify 


9? 

3u 


12 

u=0 

3u 

U=1 

3v 

u=c 

3v 

v=a 

v=b 

v=0 

(see  Figure  3.8.2) 


as  desired  outward  tangents  from  the  boundary  curves,  solving  for  the 
mixed  partial  derivatives. 

In  all  the  above  we  are  talking  about  conditions  expressed  in  terms 
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Fig.  3.8.1  Constructing  me  Surface 
from  Boundary  Curves 
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Fig.  3-8.2  Derivative  Vectors  as  Conditions 
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u=  a 
v=  t 


the  Surface 


on 


of  three-dimensional  points  or  vectors;  the  determination  of  the  homo¬ 
geneous  coordinate  representation  for  each  such  vector  is  done  by 

specifying  other  conditions,  namely,  the  point  p  . 

s 

Let  us  examine  a  particularly  elegant  way  of  determining  the 

twist  partition  of  the  surface  tensor.  In  this  analysis  we  will  use 

several  vectors  and  matrices  whose  components  are  themselves 

vectors,  such  as  the  point  00  or  the  derivative  01 

^  uv 

Let  us  require  that  a  curve  av  given  by  the  matrix  [aOal  aOyalv] 
in  endpoint  derivative  form  to  lie  across  the  surface  at  the  constant 
parameter  u=a  as  in  Figure  3.8.3.  Notice  that  as  long  as  we  pick  the 
points  aO  and  al  to  lie  respectively  on  uO  and  ul  at  the  same  value  a 
of  the  parameter  u,  the  curve  av  can  be  a  completely  arbitrary  rational 
cubic,  since  by  reparameterization  we  can  adjust  the  homogeneous  scale 
factors  at  aO  and  al  to  correspond  precisely  to  the  values  assumed  by 
the  curves  uO  and  ul  at  the  corresponding  points. 

We  thus  have 

[aO  al  aOv  aly]  =  [FQ(a)  F^a)  GQ(a)  G^a)] 


00 

01 

00 

V 

01 

V 

10 

11 

10 

V 

11 

V 

00 

u 

01 

u 

00 

uv 

01 

uv 

10 

u 

11 

u 

10 

uv 

11 

uv 

We  know  that  the  left  two  terms  [aO  al]  of  the  above  equation  will  auto- 
matically  hold  since  we  have  selected  the  points  aO,  al  to  lie  on  the 
boundaries  and  have  adjusted  the  parameterizations  properly.  The 
remaining  two  terms  give  upon  rearranging. 


[a°v  alv]  =  [F^atFjta)] 


f  00 

01 

r  oo 

01 

V 

V 

+  [G^(a)  G^(a) ] 

uv 

uv 

10 

V 

11 

v_ 

10 

L  uv 

11 

uv. 

which  further  gives 


[GQ(a)  G^a)] 


00 

uv 

01 

uv 

=  [a0v  alv]  -  [F0(a)  F^a)] 

1 

o 

o 

< 

01  1 
V 

10 

uv 

11 

uv. 

10 

L  v 

11 

vj 
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FIG.  3.8.3  INTERIOR  CURVES  AS  CONDITIONS 
ON  THE  SURFACE 
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Let  us  now  presume  another  curve  bv  =  [bO  bl  bOv  b  1  y ]  to  be 
chosen  on  the  surface,  b  #  a.  By  a  precisely  similar  argument,  we 
obtain  the  equation 


[GQ(b)  G^b)] 


poo  01  ' 

uv  uv 

10  11 

uv  uv 


=  [bOvblv]  -  [ FQ(b)  Fjtb) ] 


00  01 

v  v 

10  11 

V  V. 


Combining  the  two  equations  finally  gives 


0V 

V 


G^a)1 

"oo 

uv 

01 

uv 

aO 

V 

al 

V 

F0(a) 

Gx(b) 

10 

L  uv 

11 

uvj 

bO 

_  V 

bl 

v  J 

_F0(b) 

[00 

V 

01 

V 

10 

L  v 

11 

vj 

Since  a  or  b  are  not  0  or  1  and  are  distinct,  the  2X2  matrix 


G0(k)  Gx(a) 
LG0(b)  G^b) 


is  non- singular  and  we  can  solve  directly  for  the  twist 


partition  as 


'oo 

uv 

01 

uv 

Go(a) 

G  ^(a) 

-1 

^a0 

V 

al 

V 

"F0(a) 

FjUf 

[00 

V 

01 

V 

10 

L  uv 

11 

uvj 

LGo(b> 

G^(b) 

bO 

V 

bl 

v_ 

_F0(b) 

F^b) 

10 

V 

11 

v J 

This  result  allows  us  to  completely  characterize  a  rational  bi¬ 
cubic  surface  in  terms  of  six  essentially  arbitrary  rational  cubic 
curves  --  the  four  boundaries  and  a  pair  of  curves  passing  from  one 
boundary  to  the  opposite  boundary.  (Naturally,  the  curves  are  not  com¬ 
pletely  arbitrary  --  the  boundary  curves  must  close  and  the  two  interior 
curves  must  terminate  at  appropriate  points  on  the  boundaries.)  There 
is  remaining  one  degree  of  freedom  --  the  h^  component  mentioned 
earlier  --  which  could  be  used  to  make  the  surface  pass  through  an 
arbitrary  point,  but  at  an  unspecifiable  value  of  the  parameters  u  and  v, 
perhaps  even  outside  of  the  range  0  <  u  <  1,  0  <  v  <  1.  Solving  for  the 
value  of  h^  to  satisfy  this  condition  would  not  be  easy,  since  even 
finding  whether  a  given  surface  passes  through  a  given  point  is  a  diffi¬ 
cult  question;  numerical  approaches  would  probably  be  best. 
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3.9  PRODUCT  SURFACES  AND  SURFACES  OF  REVOLUTION 


A  convenient  way  to  specify  certain  surfaces  is  as  the  "product" 
of  two  curves;  if  the  curves  are  chosen  appropriately,  a  surface  of 
revolution  results.  Let  P(u)  and  C§(v)  be  two  curves;  then  the  i'th  com¬ 
ponent  of  a  surface  S  can  be  defined  as 


If  we  let 


S.  =  P  (u)Q.(v) 
1  11 


Pk  =  4i<u>Aik 

Qk  =  *j(v)Bjk 


where  A,  B  are  the  appropriate  matrices  for  the  basis  functions  qT,  then 
Sfc  =  <^»^(u)A^B^0j(v)  and  we  can  identify 


T...  =  A.,B.. 

ijk  lk  jk 


as  the  tensor  representing  the  surface. 

If  one  of  the  curves  ?(u)  or  ^(v)  is  planar,  say  ?(u),  we  prove 
that  all  curves  of  the  corresponding  family,  sT(u,  a),  are  planar: 

Let  X  be  the  vector  representing  the  plane  of  ?(u);  define  a 
vector  5,  component  by  component,  as 

Bi<v)  =  q^T 

Then 


S^(u,  v)B^(v)  =  P^uJQ^vlB.fv)  =  Pi(u)Ai  =  0 

where  now  the  summation  convention  holds,  and  the  curve  !?(u,  a)  lies 
in  the  plane  15(a).  Hence,  if  both  curves  ?(u),  §<v)  are  planar,  then  all 
the  constant  parameter  contours  are  plane  sections  of  the  surface. 

In  particular,  let  l5(u)  be  a  curve  in  the  x  =  y  plane  and  let  C$(v)  be 
a  portion  of  a  circle  with  center  [  0,  0,  l]  in  the  z  =  1  plane.  Then  we 
have,  in  homogeneous  coordinates. 
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pup  =  p,  p 
h  x  h  y 

(QhQx)2  +  (QhQy)2  -  «V>2  ;  QhQz  -  Qh 


as  conditions  on  the  curves.  The  surface,  in  homogeneous  coordinates, 
is 

h  -  phQh 

hx  =  PhPxQhQx 

hy  =  P,P  Q,  Q  =  P^P  Q,  Q 
J  h  y  h  y  h  x  h  y 

hz  =  PhPzQhQz  =  PhPzQh 


or,  in  regular  coordinates. 


and  we  have 


x  =  P  Q 
xx 


y  =  P  Q 
j  x^y 

z  =  P 


2  2  ,0  .2 

x  +  y  =  (pxr) 


and  also,  z  =  P  ,  giving  a  family  of  circles  about  the  z  axis.  All  plane 
z 

sections  perpendicular  to  the  z  axis  will  be  circles  and  all  plane 
sections  passing  through  the  z  axis  will  have  the  same  shape  as  the 
original  curve  P(u).  Hence,  S  is  part  of  a  surface  of  revolution  con¬ 
structed  by  rotating  the  curve  P(u)  about  the  z  axis. 


3.10  CONTINUITY  CONDITIONS 

In  constructing  an  object  from  an  assemblage  of  elementary 
surface  patches,  it  is  often  necessary  to  enforce  certain  degrees  of 
continuity  at  the  junction  between  two  patches.  In  particular,  suppose 
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S  and  T  are  two  surfaces  meeting  in  the  common  boundary  curve  B(u), 
with  the  other  boundaries  as  indicated  in  Figure  3.10.1.  Assume  the 
parameterization  in  u  and  the  homogeneous  scale  have  been  adjusted 
so  that  S(u,  1)  =  ¥(u,  0)  =  S(u),  which  is  always  possible  by  the  results 
of  subsection  2.4.  Suppose  we  wish  the  following  constraints  to  be 
satisfied: 

1.  The  Ov,  lv  boundaries  have  a  continuous  tangent  direction  at 
P  and  Q.  That  is,  the  Ov,  lv  boundary  curves  are  continuous  and  have 
a  continuous  tangent  line. 

2.  Everywhere  along  B,  the  two  surfaces  have  a  common  tangent 
plane.  That  is,  the  surfaces  are  continuous  and  have  a  continuous  unit 
normal  vector. 


It  will  turn  out  that  the  second  constraint  implies  the  first,  so  we 
will  not  explicitly  enforce  the  first  constraint. 


In  order  to  talk  about  the  tangent  planes  to  a  surface,  observe  the 
following  result. 


u 


3(ht)  _  , 


u 


t  -4-  h  t  =  h 
u 


u 


(T+ifTu) 

u  / 


^  h  “►  3 1  _► 

But  t  +  t  is  a  point  on  the  tangent  so  T^  is  a  point  on  that 

u  q-t 

tangent.  Similarly,  'f  is  a  point  on  the  tangent  -g—  .  Therefore,  a 

point  P  on  the  tangent  plane  at  "f  must  be  in  the  plane  through  'fi  'f 

and  f,  that  is,  |PTT  T  1  =  0. 

’  u  v1 

From  this  result,  the  relation  that  tangent  planes  to  surfaces  S 
and  T  coincide  can  be  expressed  as 


STT  T 
u  v 


STTT=STTT 
u  u  v1  v  u  v1 


=  0 


where  the  vectors  are  evaluated  at  the  same  value  of  the  parameter  u 
on  the  curve  B.  We  derive  this  by  noting  that  all  the  points  S’,  S  ,  S"v 
must  be  on  the  tangent  plane  to  S  and  hence  must  also  be  on  the  tangent 
plane  to  T.  By  our  initial  assumption  =  T  and  hence 

Thus  the  first  two  determinants  vanish  so  we  consider  the  third. 

In  general,  we  cannot  expect  S^  to  bear  any  relation  to  'f  or  rf 
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Fig.  3.10.1  Adjacent  Surface  Continuity  Conditions 
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so  the  vanishing  of  the  third  determinant  requires  ia  to  be  proportional 
to  f  ,  for  all  u.  But  this  relation  is  precisely  the  same  as  asking  that 
the  function  is  (u)  in  representing  a  curve  represent  the  same  curve  as 
fv(u),  regarding  the  derivatives  as  the  homogeneous  coordinates  of 
the  points  on  some  curve  in  three  dimensions.  Since  we  have  assumed 
equal  parameterization  along  the  boundary  curve  B,  this  requires 
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for  an  arbitrary  proportionality  constant  k.  We  deduce  this  by  noting 
from  the  results  of  subsection  2.4  that  a  curve  P(u)  is  the  same  as  the 
curve  Q(u)  with  the  point  P(u)  corresponding  to  the  point  Q(u)  if  and 
only  if  ?(u)  =  k^(u),  for  all  u,  for  some  constant  k. 

If  we  now  apply  an  arbitrary  reparameterization  in  u,  we  obtain 
a  more  general  set  of  sufficient  conditions  that  two  surfaces  S  and  T 
share  a  common  boundary  and  satisfy  conditions  1  and  2  as 
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where  a  =  r  =  1  will  give  equal  parameterization  along  the  boundary 
curve. 

Please  notice  that,  as  usual,  in  all  the  above  we  are  talking  about 
homogeneous  coordinate  vectors  and  that  an  equation  of  the  form 
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■01  =  k^00 
v  v 

dimensions  to  be  equal.  In  fact,  we  have 


ST^I  =  kT does  not  constrain  the  tangent  vectors  in  ordinary  three 


(S,  s)  =  (S,  )  s  +  S,  s  =  k(T,t)  =  k(T,  )  t  +  kT  t 
h  h  v  h  v  h  v  h  v  h  v 

where  s  =  [xyz  1]  and  similarly  for  t.  Now  s  =  t  but  in  general 

i-r,  ^  'f.  so  we  have 
h  h 

(S,  )  s  +  S,  ~s  =  k(T,  )  ~s  +  kT,  T 
h  v  h  h  v  h 


or 


{(S.)  -  k(T,  )  }s  =  kT,  t  -  S,  s 
h  v  h  v  >  h  v  h  v 

The  fourth  coordinate  of  this  vector  equation  gives 


and  thus 


(S,  )  -  k(T,  )  -  0 

h  v  h  v 


0  =  kT,  t  -  S,  s 
h  v  h  v 


or 


tv  kT,  Sv 
h 


Although  there  is  tangent  vector  direction  continuity  across  the 
boundary,  there  is  an  arbitrary  magnitude  discontinuity,  controlled  by 
the  constant  k. 

Notice  that  this  specific  relation  (along  with  the  corresponding 
one  for  )  forces  boundary  tangent  vector  direction  continuity  at  P 
and  Q,  satisfying  condition  1 
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SECTION  IV 


IMPLEMENTATION 


4.1  EXAMPLES 

Figures  4.1.1  through  4.1.5  have  been  included  to  illustrate  the 
appearance  of  the  surfaces.  The  figures  are  all  true  perspective  pro¬ 
jections  of  the  surfaces  from  varying  vantage  points.  They  should  be 
viewed  from  such  a  distance  as  to  give  a  45°  field  of  view.  They  were 
computed  upon  a  DEC  PDP-1  with  DEC  340  scopes  used  for  display. 

The  PDP-1  is  a  rather  ancient  machine  by  modern  standards,  having 
an  18-bit  word  length,  5-qsec  cycle  time,  25-qsec  multiply,  40-qsec 
divide  and  no  floating  point. 

The  time  to  compute  a  new  view  of  a  surface  is  1  to  2  seconds 
for  these  figures,  depending  on  the  number  of  curves  used.  This  is 
about  1-2  ms.  per  chord.  Each  curve  is  approximated  by  32  chords. 

The  tensors  for  the  cylinder,  sphere  and  toroid  were  manually 

computed  as  surfaces  of  revolution.  The  hyperbolic  paraboloid  or 

2  2 

saddle  surface  is  simply  the  equation  z  =  x  -  y  expressed  parametri¬ 
cally  in  the  tensor.  The  surface  shown  in  Figure  4.1.5  was  obtained  by 
making  a  20%  random  variation  in  each  of  the  components  of  the  tensor 
for  the  cylinder  of  Figure  4.1.2. 

4.2  EXPERIMENTAL  IMPLEMENTATION 

Faced  with  the  speed  and  capacity  restrictions  of  the  available 
interactive  computing  facilities  (the  PDP-1  graphics  research 
laboratory),  we  decided  that  the  best  way  to  illustrate  the  capabilities 
of  this  material  would  be  the  experimental  production  of  a  short  motion- 
picture  film  showing  representative  examples  of  the  surfaces  moving 
about  to  give  a  good  impression  of  their  three-dimensional  structure 
and  appearance.  This  section  comments  upon  the  results  of  that  experi¬ 
ment,  with  the  material  in  Appendix  II  proposing  in  relation  to  this 
section  what  we  would  like  to  have  tried,  given  a  sufficiently  more 
powerful  interactive  environment. 
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Fig.  4.  1. 2  Cylinder 
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Fig.  4.1-4  Sphere 
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4.1.5 


Modified  Cylinder 


Two  concerns  were  foremost  in  the  design  and  coding  of  the 
routines  for  making  the  film.  The  first  was  to  incorporate  the  theory 
and  algorithms  into  a  moderately  flexible  system  in  which  as  much 
responsive  interaction  as  possible  could  be  used  to  generate  an  aes¬ 
thetic  sequence  of  continuously  moving  views  of  the  surfaces.  Oppor¬ 
tunity  for  experimentation  with  the  composition  and  form  was  to  be 
provided  on  the  computer  long  before  the  sequence  was  committed  to 
film.  The  second  concern  was  to  make  the  resulting  program  run  as 
fast  as  possible  during  the  crucial  phase  of  generating  a  visual  image 
of  the  surface.  This  was  accomplished  by  using  few  subroutines  in  the 
critical  areas,  writing  in  machine  language,  and  by  unwinding  loops 
into  straight-line  code.  Subscripted  arrays  had  to  be  removed  since 
the  PDP-1  has  no  index  registers  and  must  compute  variable  subscripts 
with  the  standard  arithmetic  instructions,  using  indirect  addressing  for 
references.  As  a  compromise  forced  by  the  relatively  slow  speed  of 
the  computer  and  displays,  magnetic  tape  was  used  to  store  some  inter¬ 
mediate  information. 


The  surfaces  considered  in  the  system  are  time-varying  surfaces 
of  the  form 


S 


ijk 


+  (l-a)S’’k 


where  the  "mixing"  parameter  a,  0  <  a  ^  1,  is  continuously  adjustable 
from  a  potentiometer.  Thus  we  can  have  a  surface  which  changes  from 
a  plane  to  a  hemisphere  in  the  process  of  continuous  deformation.  Very 
interesting  intermediate  surfaces  are  formed  in  such  cases  as  trans¬ 
forming  a  hemisphere  into  a  quarter -torus. 


The  operation  of  the  system  is  shown  in  Figure  4.2.1  and 
described  as  follows.  The  tensors  for  a  set  of  possible  surface  patches 
had  been  computed  by  hand  very  early  in  this  research  and  were  now 
made  known  to  the  system.  Compound  surfaces  formed  from  pairs  of 
patches  are  defined  as  the  objects  of  interest.  Several  such  compound 
surfaces  can  be  displayed  at  once,  allowing  a  sphere  to  be  formed  from 
two  hemispheres  or  reference  planes  to  be  added.  The  compound 
surface  also  contains  a  specification  of  the  appearance  it  is  to  take  on 
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FIG.  4.2.1  FILM  MAKING  SYSTEM 
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the  display,  that  is,  the  number  of  parametric  curves  and  the  number  of 
points  on  the  curves.  For  initial  composition,  a  simple  visual  repre¬ 
sentation  --  say  three  curves  for  each  parameter,  five  points  for  each 
curve  --  is  displayed.  Potentiometers  are  used  to  control  the  rates  of 
rotation  of  the  surface  about  the  three  different  axes  and  to  control  the 
rate  at  which  the  compound  surfaces  metamorphose  from  one  case  to 
another.  Translation  along  the  three  axes  is  possible  by  a  set  of  switch 
controls,  but  this  type  of  motion  did  not  add  much  to  the  film  sequence. 
During  this  continuous  motion  the  successive  positions  in  time  -- 
represented  by  a  4  X  4  homogeneous  matrix  --  and  the  mixing  parameter 
are  recorded  on  tape.  These  computations  for  a  simple  representation 
can  be  performed  at  a  rate  of  about  five  to  ten  frames  per  second,  suf¬ 
ficiently  close  to  real-time  to  make  the  potentiometer  controls  adequate 
for  composing  the  film  sequence.  The  visual  representation  can  then  be 
changed  to  a  more  suitable  detail,  say  17  curves  and  33  points  per  curve. 
The  position  tape  is  used  as  data  to  generate  the  actual  frames  of  the 
film  for  this  detailed  representation,  which  are  recorded  on  another  tape. 
This  computation  is  considerably  slower,  on  the  order  of  one  to  five 
seconds  per  frame,  depending  on  the  number  of  surfaces  displayed  at 
once  and  the  fineness  of  the  visual  representation.  The  perspective  view¬ 
point  is  chosen  at  the  same  time  as  the  visual  representation  and  can  be 
varied  for  different  instances  of  the  same  sequence  of  positions.  The 
final  visual  tape  can  be  read  at  the  maximum  speed  allowed  by  the  tape 
drives  and  the  display  processor  --  about  ten  frames  per  second  --  to 
let  the  detailed  sequence  be  viewed  (but  not  modified)  in  almost  real¬ 
time  without  recording  on  film.  The  images  on  the  tape  are  then  photo¬ 
graphed  a  frame  at  a  time  on  a  high-quality  35 -mm  animation  camera. 

As  is  probably  true  in  many  interactive  systems,  the  need  to 
achieve  fast  response  forced  much  time  to  be  spent  on  the  rather  unil¬ 
luminating  problems  of  writing  efficient  machine  code  rather  than  on 
the  more  important  areas  of  human-engineering  and  overall  aesthetic 
design.  These  problems  were  of  course  aggravated  by  the  nature  of  the 
objects  being  computed  --  realistic  representations  of  three-dimensional 
surfaces.  It  is  clear  that  computer  graphics  is  one  major  area  of 
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computational  endeavor  in  which  even  the  simplest  of  problems  can  be 
expected  to  reach  the  bounds  of  the  available  computing  resources.  It 
is  also  true  that  relatively  minor  improvements  can  be  expected  to  have 
major  effects  on  the  practicality  of  a  proposed  procedure.  In  the  case 
of  these  bi-cubic  rational  surfaces,  either  a  special  purpose  curve 
generator  as  described  above  in  subsection  2.2  or  a  conventional 
general  purpose  processor  with  a  few  hundred  words  of  fast  memory  to 
hold  the  subroutines  for  computing  a  curve  (including  the  homogeneous 
coordinate  division  and  generation  of  display  processor  code)  would 
enable  a  useful  real-time  system  for  the  drawing  of  space  forms  using 
the  representation  discussed  here.  Even  more  simply,  a  highly  opti¬ 
mizing  compiler  for  an  algebraic  language  with  a  library  of  standard 
matrix  routines  would  have  reduced  the  programming  time  for  a  substan¬ 
tial  part  of  the  movie  system  to  a  matter  of  hours,  rather  than  the 
month  or  so  it  took.  If  the  compiler  had  facilities  such  as  the  proposed 
I.B.M.  PL/l  compile-time  macros,  it  would  be  possible  to  expand  the 
inner  loops  and  achieve  as  fast  running  a  program  as  could  be  expected, 
even  on  a  computer  like  the  PDP-1  with  a  simple  fixed-point  instruction 
set. 

Note:  Copies  of  a  twelve -minute,  16-mm,  black  and  white  silent 
film  made  by  this  process  may  be  obtained  by  purchase  order  to 

Film  Service  Laboratory 
62  Berkely  Street 
Boston,  Massachusetts  02116 
Telephone:  (61 7)-542-1238 

In  ordering  please  request  a  reversal  print  of  "Rational  Bi-Cubic 
Surfaces"  by  Theodore  M.  P.  Lee,  Harvard  University.  The  price  will 
be  about  $40.  The  film  includes  examples  of  all  the  surfaces  in 
Appendix  I,  each  surface  slowly  changing  into  another  surface.  The  film 
is  to  be  projected  at  the  standard  sound  speed  of  24  frames  per  second. 
(Because  of  idiosyncracies  in  the  35  mm  to  16  mm  reduction  process,  a 
non-reversal  print  gives  an  annoying  non-symmetric  white  border  around 
the  frame  --  the  reversal  print  will  have  black  lines  on  a  white  background 
and  will  have  higher  contrast  and  definition  than  the  non-reversal  print.) 
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SECTION  V 


PROBLEMS  FOR  FURTHER  RESEARCH 

5.1  PROBLEMS  FOR  FURTHER  RESEARCH 

The  following  paragraphs  outline  a  few  of  the  important  questions 
open  for  future  research.  I  cover  here  only  the  more  mathematical 
questions;  clearly  the  problems  associated  with  organizing  these  results 
into  an  effective,  human-engineered  system  for  computer  aided  design 
must  be  answered  at  some  time. 

Since  a  boundary  curve  can  be  degenerate,  patches  with  only  two 
or  three  sides  are  possible.  Such  patches  by  virtue  of  the  degenerate 
boundary  curve  must  be  treated  differently  as  far  as  tangent  plane  conti¬ 
nuity  is  concerned.  Exactly  what  constraints  can  be  imposed  and  what 
types  of  surfaces  with  degenerate  curves  exist  are  open  questions. 

Another  question  of  continuity  constraint  is  the  general  problem 
of  two  overlapping  surfaces,  say  an  airplane  wing  and  fuselage.  Here  we 
are  attempting  to  ask  for  surface  continuity  in  the  middle  of  a  patch,  not 
on  the  boundary.  Of  course,  we  could  use  contours  in  the  interior  of  a 
patch  to  define  a  new  subpatch  at  the  boundaries  of  which  continuity  is  to 
be  imposed.  But  this  may  not  always  be  possible,  nor  may  it  be  the 
most  general  solution. 

Further  results  are  needed  in  specifying  surfaces  from  external 
data.  Personally,  I  prefer  the  surface-molding  approach  in  which  on¬ 
line  feedback  is  continuously  present  as  contrasted  to  the  point  surface¬ 
fitting  approach  currently  used  by  most  surface  design  systems.  How¬ 
ever,  there  are  many  occasions  in  which  it  would  be  expedient  to  use 
such  easily  displayed  surfaces  as  these  to  approximate  a  surface  de¬ 
scribed  in  some  other  form,  perhaps  from  measured  coordinates 

Even  using  these  bi-cub ic  rational  surfaces  as  molded  surfaces, 
there  are  not  at  present  sufficient  techniques  for  manipulating  the  dis¬ 
play.  For  instance,  efficient  methods  of  performing  small  variations  in 
the  surface  are  obviously  needed. 
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I  will  do  no  more  than  mention  the  question  of  discovering  hidden 
lines  in  a  figure  composed  of  these  surfaces.  The  problem  for  this 
geometry  is  obviously  very  difficult  and  defies  a  simple  analytical 
solution.  Even  the  application  of  Warnock's  algorithm  [23],  effective  as 
it  is,  may  not  be  the  answer,  since  his  techniques  have  been  found  very 
susceptible  to  Mach  band  distortion  on  curved  surfaces  approximated  by 
planar  sections.  Such  distortion  gives  the  surface  a  "fluted"  appearance 
at  the  joints  of  the  planes  caused  by  a  psychological  enhancement  of  the 
intensity  discontinuity. 
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APPENDIX  I 


EQUATIONS  OF  TYPICAL  OBJECTS 

We  list  here  the  equations  of  six  typical  surfaces.  They  were  used 
for  the  photographs  in  Section  IV  and  as  the  basic  surfaces  used  in  the 
movie  described  in  subsection  4.2.  The  sphere,  cone,  cylinder  and 
bottle  describe  one -half  of  the  surface  by  that  name;  the  torus  describes 
one-quarter  of  a  complete  torus.  The  bottle  surface  was  formed  as  a 
surface  of  revolution  whose  generator  is  the  function  F^;  it  looks  vaguely 
like  a  bottle  neck  or  wine  glass  stem. 

1.  Cylinder 
hx  =  2v  -  1 

hy  =  -2v^  +  2v 
hz  =  u(  2v^-2v+l) 
h  =  2v^  -  2v  +  1 

2.  Cone 

hx  =  (2u-l)v 
hy  =  (-2u^+2u)v 
hz  =  v(2u^-2u+l) 
h  =  2u2  -  2u  +  1 

3.  Saddle 

hx  =  2u  -  1 
hy  =  2v  -  1 

hz  =  4u^  -  4u  -  4v^  +  4v 
h  =  1 
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4.  Sphere 


hx  =  (2u-l)(-2v2+2v) 
hy  =  (-2u2+2u)(-2v2+2v) 
hz  =  (2u2-2u+1)(2v-1) 
h  =  (2u2-2u+l)(2v2-2v+l) 

5.  Torus 

hx  =  (36u2  -  36u+  9)(6v  -  3) 
hy  =  (36u2-36u+9)(18v2-  18v+4) 
hz  =  (6u  -  3)(18v2  -  18v+  5) 
h  =  (18u2  -  18u+  5)(18v2  -  18v+  5) 

6.  Bottle 

hx  -  (2u  -  l)(2v2  -  3v“  +  1) 
hy  =  -(-2u2  +  2u)(2v2  -  3v2  +  1) 
hz  =  2(2u2  -  2u+  l)(v  -  ) 

h  =  2(2u2  -  2u+  1) 
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APPENDIX  II 


PROPOSAL  FOR  A  DATA  STRUCTURE 

The  mathematics  of  Sections  II  and  III  provides  a  set  of  tools  to 
be  used  in  the  specification  and  display  of  curves  and  surfaces.  These 
results  provide  the  muscles  for  a  computer-aided  design  system  but 
neither  the  skeletal  framework  holding  it  together  nor  the  nervous 
system  enabling  and  controlling  action;  neither  do  they  indicate  the 
ways  in  which  the  system  is  to  relate  to  the  world,  nor  the  means  and 
direction  of  its  evolution  and  growth.  In  this  appendix  we  propose  a 
data  base  structured  so  as  to  serve  as  a  framework  into  which  the 
mathematics  will  naturally  fit  and  onto  which  it  will  operate.  This  data 
base  is  to  contain  explicit  models  of  the  connections  between  the  con¬ 
cepts  involved  in  talking  about  surfaces  and  is  designed  to  expand  quali¬ 
tatively  as  these  concepts  are  expanded. 

This  material  was  developed  in  response  to  the  difficulties 
encountered  in  programming  the  examples  described  in  Section  IV.  It 
became  clear  that  a  significantly  useful  implementation  of  the  curve 
and  surface  theory  would  require  the  development  of  routines  that 
really  "knew"  about  surfaces  rather  than  ones  that  were  only  capable 
of  blindly  manipulating  the  tensors  in  the  fashion  of  the  present 
programs.  The  primary  evidence  for  this  conclusion  was  the  aware¬ 
ness  that  throughout  the  development  and  use  of  the  programs, 
experience  would  suggest  continual  changes  in  their  mathematical  and 
functional  appearance.  These  changes  need  to  be  incorporated  into  the 
system  in  as  painless  a  manner  as  possible.  A  tentative  solution  to  the 
problem  was  suggested  by  the  extreme  generality  and  flexibility  of 
Feldman's  associative  data  structure  [8]  Since  none  of  the  material  in 
this  section  has  been  tested,  it  is  presented  only  as  groundwork  that 
should  not  have  to  be  retraced  in  the  future. 

The  intent  of  this  material  is  to  demonstrate  the  necessity, 
practicality  and  utility  of  a  highly  general  data  structure  in  the  particu¬ 
lar  context  of  computer-aided  design  of  space  forms.  We  will  pursue 
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the  details  only  so  far  as  is  necessary  to  make  these  points  and  to 
illustrate  the  possibilities  inherent  in  the  homogeneous  tensor  repre¬ 
sentation  chosen  for  curves  and  surfaces.  In  this  presentation  there 
will  be  many  design  questions  --  such  as  the  manner  in  which  matrices 
are  to  be  represented  --  whose  resolution  will  depend  upon  experiment. 
In  these  areas  I  am  being  purposely  vague  in  order  not  to  make  an 
unfounded  commitment  to  some  particular  implementation.  I  make  no 
claim  of  originality  for  the  immediately  following  general  comments  on 
data  structures,  although  I  can  cite  no  particular  source  for  them; 
their  sole  purpose  is  to  set  the  stage  for  the  case  under  discussion. 

Philosophically  speaking,  the  purpose  and  appropriateness  of  a 
data  structure  are  operationally  determined  by  the  uses  to  which  it  is 
put,  that  is,  by  the  questions  the  data  structure  is  expected  to  help 
answer  and  what  those  answers  are  expected  to  accomplish.  We  think 
of  the  data  structure  as  being  used  to  form  a  model  of  some  possibly 
changing  set  of  concepts;  instead  of  asking  questions  about  the  concepts 
directly,  we  ask  questions  about  the  model  and  expect  that  the  answers 
can  be  interpreted  as  answers  to  equivalent  questions  asked  about  the 
concepts. 

The  choice  of  a  data  structure  thus  implies  not  only  a  particular 
modelling  of  the  concepts  but  also  a  modelling  of  the  questions:  to 
answer  a  question  put  to  the  set  of  concepts,  we  map  the  question  to  its 
appropriate  model,  ask  that  question  or  questions  of  the  data  base, 
built  in  accordance  with  the  chosen  data  structure,  and  map  the  answer 
back  to  an  answer  of  the  original  question.  If  we  choose  a  simple  data 
structure  to  model  a  complex  set  of  concepts,  a  simple  question  asked 
in  the  language  of  the  set  of  concepts  will  be  mapped  into  a  complex 
question  asked  in  the  language  of  the  data  structure.  If  either  the 
resulting  question  or  the  mapping  determining  it  --  presumably  some 
program  --is  too  complex,  it  is  obvious  that  the  data  structure  is 
unsuited  to  finding  an  answer  to  the  original  question;  in  fact,  it  may 
be  impossible  to  do  so.  In  choosing  a  data  structure,  if  we  do  not  know 
in  advance  all  the  questions  to  be  asked,  we  at  least  ought  to  be  able  to 
convince  ourselves  that  the  types  of  questions  likely  to  arise  will  not  be 
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too  difficult,  the  difficulty  being  measured  in  terms  of  the  mappings 
necessary  to  apply  the  question  to  the  structured  data  base. 

Now  let  us  talk  about  curves  and  surfaces.  The  claimed  result  of 
this  dissertation  is  that  the  homogeneous  tensor  representation  of 
surfaces  as  rational  bi-cubics  is  desirable.  We  thus  assume  that  part 
of  the  representation  of  a  surface  will  be  such  a  tensor.  Immediately 
we  notice  that  this  tensor  is  not  unique,  since  there  are  several  differ¬ 
ent  such  representations  suitable  for  different  purposes: 

1.  The  polynomial  representation  is  good  as  an  abstract  proto¬ 
type  of  the  surface. 

2.  The  endpoint  derivative  representation  is  good  for  specifying 
the  curve  from  geometrical  considerations. 

3.  Any  of  several  difference  equation  representations  is  good  for 
generating  a  visual  representation  of  the  surface,  the  particular 
difference  equation  used  dependent  on  how  many  and  what  type  of  curves 
we  wish  to  display. 

4.  The  code  driving  the  display  is  itself  a  representation,  making 
the  surface  visible. 

In  addition  to  these  explicit  representations  applicable  to  the 
general  bi-cubic  rational  surface,  we  have  the  representation  of  a 
product  surface  in  terms  of  two  curves,  each  curve  in  turn  represented 
as  a  matrix.  And  then  we  have  a  surface  of  revolution  represented  as  a 
product  surface  in  which  one  of  the  curves  is  an  arc  of  a  circle  and  the 
other  an  arbitrary  plane  curve.  Or,  by  the  results  of  subsection  3.8,  we 
can  represent  a  surface  by  six  curves  and  a  few  additional  parameters. 
In  the  course  of  working  with  a  particular  surface,  all  of  these  different 
representations  might  be  useful  and  might  exist  at  the  same  time.  We 
can  also  conceive  of  additional  representations  or  characterizations  as 
yet  undeveloped  in  detail,  such  as  requiring  that  a  surface  be  a  torus,  a 
sphere,  the  hull  of  a  three-masted  schooner,  or  the  right  tail-fin  of  an 
Edsel  station  wagon,  assuming  all  such  concepts  can  be  formulated  to  a 
sufficient  degree  of  exactness.  We  will  not  succeed  in  building  an 
interesting  surface  manipulation  system  unless  it  is  capable  of  freely 
wandering  about  this  proliferation  of  representations,  choosing  the 
appropriate  one  for  the  task  at  hand. 
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If  we  grant  that  any  or  all  of  these  various  representations  of  a 
surface  need  to  be  modelled,  we  see  two  types  of  questions  to  be 
answerable  by  use  of  the  data  base.  Prototypes  of  these  questions  are: 

1.  "What  is  the  display  code  to  be  used  to  draw  this  surface?" 

and 

2.  "I  find  that  there  is  no  display  code  for  this  surface;  how  do 
I  generate  it  in  order  that  I  may  answer  question  1?" 

The  type  of  answer  to  the  first  question  is  clear;  to  the  second 
question  we  would  expect  an  answer  of  the  form:  "if  you  have  a  finite 
difference  representation  of  the  surface,  you  can  use  subroutine 
GENERATESURFACE  to  compute  the  display  code."  Notice  that  this 
answer  implies  the  additional  question  --  "What  is  the  finite  difference 
representation  of  the  surface?"  --  which  may  lead  to  further  questions 
about  finite  difference  representations.  In  the  use  of  most  common 
data  structures,  it  is  only  questions  of  type  (1)  that  are  able  to  be 
answered  directly;  questions  of  type  (2)  are  treated  implicitly  by  the 
programs  that  massage  the  data  base.  Changing  these  programs  in 
order  to  expand  this  set  of  answerable  questions  is  usually  very  painful. 

Another  way  to  characterize  the  two  questions  is  to  notice  that  in 
the  first,  explicit  information  answers  the  question;  in  the  second, 
explicit  information  answering  the  question  about  finite  difference 
representations  yields  the  implicit  answer  to  the  question  about  display 
code  --  that  is,  implicit  in  the  fact  that  the  object  has  a  finite  dif¬ 
ference  representation  is  the  fact  that  display  code  for  that  object  can 
be  generated  by  a  particular  subroutine. 

If  we  permit  the  data  structure  to  allow  modelling  of  the  type  of 
information  implied  in  question  (2),  we  can  provide  an  answer  to  newly 
arising  questions  in  terms  of  previously  known  questions.  For  instance, 
if  we  define  a  sphere  as  a  particular  kind  of  surface,  describing  the 
way  to  obtain  a  tensor  representation  of  the  surface  from  some  particu¬ 
lar  simple  representation  of  the  sphere,  we  can  answer  any  question 
about  the  sphere  if  the  same  question  can  be  answered  about  a  general 
surface.  Continuing,  if  we  now  decide  that  it  is  useful  to  talk  about  the 
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boundary  curves  of  a  surface,  specifying  how  that  is  to  be  done  in  terms 
of  the  tensor  representation  of  a  general  surface,  we  can  ask  for  the 
boundary  curves  of  a  sphere,  even  though  the  sphere  is  initially  repre¬ 
sented  only  by  its  center  and  radius,  having  no  explicitly  represented 
boundary  curves. 

Abstracted  from  the  above,  the  types  of  information  that  need 
direct  representation  are  as  follows: 

1.  Explicit  pieces  of  numeric  information,  such  as  the  contents 
of  a  tensor  or  matrix. 

2.  Objects  in  the  abstract,  such  as  a  particular  surface  or  curve, 
apart  from  any  specific  numeric  characterization. 

3.  Relations  between  objects,  such  as  "this  tensor  is  the  poly¬ 
nomial  representation  of  that  surface,  "  or  "this  curve  is  a  boundary 
curve  of  that  surface." 

4.  Generic  information,  that  is,  the  answers  to  questions  of 
type  (2).  We  will  suggest  a  way  in  which  for  this  problem  area  this  type 
of  information  can  be  modelled  in  terms  of  the  first  three. 

Categories  (1)  through  (3)  are  to  be  structured  to  allow  direct 
answers  to  questions  of  type  (1)  --  given  an  object,  any  of  the  pertinent 
properties  of  the  object  are  intended  to  be  retrieved  by  examining  the 
appropriate  relations  holding  between  the  object  and  the  desired  property. 
(This  will  be  made  clearer  shortly.)  Any  of  a  number  of  schemes  for 
representing  relations  could  be  chosen  so  as  to  permit  the  direct  access 
to  this  kind  of  information.  We  will  briefly  mention  later  one  such 
scheme,  but  our  chief  concern  is  to  point  out  in  this  context  the  utility 
of  the  general  concept  of  representing  information  in  terms  of  relations. 

To  amplify  category  (3)  recall  [9,  p.  49]  the  formal  definition  of 
a  relation  as  a  set  of  ordered  pairs  drawn  from  some  universe,  such  as 
=  =  {(  x,  y)|  x  2  yj .  An  instance  of  a  relation  is  just  one  of  the  ordered 
pairs  --  (  2,  3)  --  which  is  better  expressed  as  the  ordered  triple  (2^3) 
in  order  to  emphasize  the  relation  of  interest.  Typically,  we  will  be  con¬ 
cerned  with  mathematical  relations  such  as  set  membership  or  with 
practical  relations  such  as  the  one  holding  between  a  surface  and  any  of 
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its  several  representations.  Both  an  entire  relation  and  a  specific 
instance  of  a  relation  must  be  considered  as  objects  in  themselves, 
which  can  participate  in  relations  with  other  objects.  In  this  way  it 
would  be  possible  to  indicate,  for  instance,  that  the  relation  ^  is 
transitive  or  that  the  particular  statement  "The  finite  difference 
representation  of  the  surface  A  is  B"  was  derived  by  using  generic 
information  about  finite  difference  representations 

If  we  thus  talk  explicitly  about  relations  and  represent  them 
explicitly  in  the  data  base,  it  becomes  possible  to  introduce  new 
concepts  --  expressed  as  relations  --  and  to  relate  those  new  con¬ 
cepts  to  existing  concepts  in  an  operationally  useful  manner.  This  is, 
of  course,  not  a  new  idea;  I  present  it  here  with  the  claim  that  the 
problems  here  discussed  are  well-suited  to  solution  by  these  means. 

At  this  point  it  seems  best  to  illustrate  what  I  mean  by  an 
extended  set  of  examples.  These  examples  are  contained  in  the  graphs 
in  Figures  A. 2.1  through  A. 2. 12  and  will  be  discussed  in  detail  in  the 
next  several  pages.  The  graphs  are  directed,  with  labels  (possibly 
empty  when  unessential)  on  both  nodes  and  arcs.  The  graphs  indicate 
some  portion  of  the  information  to  be  represented  in  the  data  base  and 
thus  consist  of  a  number  of  relations  between  objects,  a  relation  being 
defined  by  the  explicit  representation  of  its  instances.  The  label  on  an 
arc  or  node  is  intended  to  be  a  unique  name  for  the  object  represented 
by  the  arc  or  node  The  instance  (ARB)  of  a  relation  R  is  indicated  by 
making  an  arc  labelled  R  from  a  node  labelled  A  to  a  node  labelled  B. 
The  fact  that  an  arc  and  a  node  may  have  the  same  name  is  a  result  of 
allowing  relations  to  be  treated  as  objects;  in  order  to  treat  an  instance 
of  a  relation  as  an  object,  we  will  circle  the  label  appearing  on  the  arc 
representing  the  instance  and  treat  that  circle  as  if  it  were  a  node  in 
the  graph. 

As  our  first  example,  let  us  look  at  the  information  intended  to 
be  retrieved  by  the  prototype  questions  (1)  and  (2)  above.  If  we  denote 
the  particular  surface  concerned  by  S,  the  first  question  translates  to 
the  request  "Find  x  such  that  the  relation  (S  DISPLAY-CODE  x)  holds." 
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The  answer  expected  from  this  question  would  be  found  in 
(S  DISPLAY-CODE  D)  as  the  object  D.  The  second  question  can  be 
rephrased,  not  quite  so  formally,  as  "is  there  a  statement  in  the  data 
base  which  could  be  used  to  deduce  the  existence  of  a  subroutine  to  com¬ 
pute  display  code  for  an  object,  and  if  so,  can  that  deduction  be  made  for 
the  particular  object  S?"  The  answer  to  this  question  would  be  supplied 
by  a  representation  of  the  statement  "For  all  surfaces  x  the  existence 
of  a  finite  difference  representation  y  implies  the  existence  of  display 
code  z  whose  value  is  computed  from  y  by  the  subroutine 
GENERATEDISPLAY." 

Figures  A. 2.1  and  A. 2. 2  are  two  models  for  this  statement;  an 
exegesis  of  the  models  follows.  To  begin,  the  statement  can  be  written 
formally  as  the  conditional  proposition 

(V x,  y)  (x  e  SURFACES)  (x  FIN.  DIF.  REP.  y)  =4 

(3t,  z)  (x  DISPLAY-CODE  z)  (t  ARG.  y)  (t  SUBR.  GENERATEDISPLAY)  (t-*z) 

where  the  object  t  represents  the  application  of  the  subroutine,  the 
relation  ARG.  describes  y  as  an  argument  of  the  subroutine  and  the 
relation  SUBR.  gives  the  name  of  the  subroutine  to  be  called;  the  instance 
(t-*>z)  says  that  this  subroutine  call  will  give  the  value  of  z. 

Figure  A. 2.1  represents  the  statement  by  labelling  the  four  nodes 
corresponding  to  the  quantified  variables  x,  y,  z,  t  with  the  appropriate 
quantification  symbol;  the  labels  for  the  variables  are  omitted  since 
they  are  irrelevant  and  only  serve  to  distinguish  the  separate  variables, 
which  is  done  here  by  using  separate  nodes.  The  two  relations  serving 
as  the  antecedents  of  the  conditional  statement  are  drawn  with  heavy 
lines;  the  four  consequents  are  drawn  with  dotted  lines.  In  general,  we 
read  a  diagram  such  as  this  by  saying  "if  all  the  relations  indicated  by 
heavy  lines  are  true,  for  appropriate  substitutions  for  the  quantified 
variables,  then  we  can  deduce  that  all  the  relations  indicated  by  dotted 
lines  are  true." 

To  actually  represent  this  statement  strictly  in  terms  of  objects 
and  relations,  we  would  have  to  use  a  scheme  like  that  in  Figure  A. 2. 2. 
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FIG.  A.  2.1  REPRESENTATION  OF  DISPLAY  CODE 

DEFINITION  (SIMPLIFIED) 
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SURFACES 


FIG.  A.  2.2  REPRESENTATION  OF  DISPLAY  CODE  DEFINITION 

(DETAILED) 


The  node  P  groups  the  set  of  relations  and  variables  in  the  statement 
together;  the  part  played  by  each  object  or  instance  of  a  relation  in  the 
statement  is  indicated  by  relating  that  object  or  instance  to  P.  Quanti¬ 
fied  variables  are  related  by  the  quantification  symbol,  consequents  by 
the  relation  CONS,  and  antecedents  by  ANT. 

We  would  apply  this  representation  to  answering  question  (2)  by 
translating  question  (2)  into  a  sequence  of  questions  of  type  (1),  namely, 

(a)  "Find  existentially  quantified  variables  x,  t  and  a  universally 
quantified  variable  z  such  that  the  instances  (x  DISPLAY-CODE  z)  and 
(t  -*■  z)  are  consequents  of  some  statement  P." 

(b)  "For  every  antecedent  of  the  statement  P,  verify  that  the 
antecedent  relations  hold  when  the  universally  quantified  variables  are 
replaced  by  the  appropriate  particular  objects  in  the  data  base,  such  as 
by  substituting  the  particular  surface  S  for  the  variable  x.  If  all  the 
antecedent  relations  hold,  create  new  objects  for  the  existentially  quanti¬ 
fied  variables  and  enter  them  into  the  data  base  via  the  consequent 
relations." 

We  will  not  explore  this  deductive  process  much  further;  the 
literature  on  deductive  question-answering  systems  [3,  25,  and  works 
there  cited]  seems  to  indicate  the  feasibility  of  this  approach  and 
covers  in  much  more  detail  than  I  care  to  know  about  at  the  present  the 
precise  algorithms  to  follow.  I  claim  without  other  than  empirical 
justification  that  the  form  of  the  statement  represented  by  these  figures 
is  sufficient  to  answer  the  questions  of  type  (2)  that  arise  in  this  appli¬ 
cation.  The  subsequent  examples  are  all  of  this  type  --  conditional 
expressions  indicating  what  relations  can  be  deduced  from  what  relations 
--  and  in  a  few  pages  cover  a  good  cross  section  of  the  concepts 
developed  in  the  previous  sections. 

(As  an  aside,  the  fact  that  the  types  of  deduction  to  be  used  in 
accessing  the  proposed  data  base  can  be  limited  to  this  relatively 
simple  form  ought  to  point  heavily  toward  the  likelihood  of  a  practical 
system  based  on  these  ideas.) 
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Before  illustrating  the  representation  of  some  specific  concepts 
about  surfaces,  we  introduce  representations  of  three  mathematical 
concepts  that  will  simplify  the  subsequent  representation  of  the  surface 
material.  The  concepts  to  be  discussed  are 

1.  set  membership 

2.  change  of  basis  transformations,  and 

3.  matrix  multiplication. 

Representing  these  concepts  directly  in  the  data  base  rather  than  merely 
using  them  in  programs  serves  both  the  purpose  of  forcing  the  concepts 
out  in  the  open  where  they  will  be  subject  to  careful  definition  and  the 
purpose  of  enabling  their  application  to  a  wide  variety  of  circumstances. 

In  Figure  A. 2. 3  is  represented  the  statement  that  a  member  of  a 
set  is  a  member  of  any  including  set;  that  is,  for  any  objects  x,  y  and 
z,  given  (x  e  y)  and  (y  C  z),  we  can  deduce  (x  e  z).  Having  introduced 
this  concept  we  no  longer  need  to  completely  characterize  objects  which 
are  defined  as  belonging  to  subclasses  of  existing  classes;  for  instance, 
having  defined  the  class  of  spheres  to  be  a  subset  of  the  class  of  surfaces, 
any  deductions  applicable  to  an  object  S  bearing  the  explicit  relation 
(S  e  SURFACES)  will  be  applicable  to  an  object  T  bearing  the  relation 
(T  e  SPHERES),  assuming  any  other  antecedent  conditions  also  hold. 

In  Figure  A. 2.4  is  defined  once  and  for  all  the  inclusion  relation¬ 
ships  among  the  classes  of  interest.  Having  done  so  we  do  not  need  to 
define  all  the  properties  of  each  class  of  objects  and  can  concentrate  on 
the  essential  properties  peculiar  to  that  class. 

The  essential  differences  between  the  various  tensor  represen¬ 
tations  of  surfaces  can  be  phrased  in  terms  of  their  different  functional 
bases.  As  mathematicians  we  know  how  to  transform  a  tensor  repre¬ 
sented  in  one  basis  to  a  tensor  represented  in  a  different  basis;  namely, 
if  there  is  a  transformation  between  two  bases,  a  representation  of  an 
object  in  the  range  basis  is  the  matrix  product  of  the  representation  in 
the  domain  basis  and  the  transformation  matrix.  Figure  A. 2. 5  is  a 
diagrammatic  version  of  this  statement  --  its  precise  interpretation 
would  be  that  the  indicated  matrix  product  relates  the  two  representations. 
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FIG-  A.  2.3  SET  MEMBERSHIP 


CYLINDERS 


FIG.  A.  2.4  HIERARCHY  OF  OBJECTS 
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FIG.  A.2.5  TRANSFORMATION  OF 
BASIS 
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Hereafter,  to  enable  the  conversion  of  a  representation  of  a  surface 
from  one  basis  to  another,  we  only  need  define  the  particular  transfor¬ 
mation  matrix  relating  the  bases. 

In  the  above  paragraph  we  have  defined  the  change  of  basis  trans¬ 
formation  in  terms  of  a  matrix  product,  representing  the  matrix 
product  symbolically.  In  Figure  A. 2. 6  we  represent  the  additional  fact 
that  some  subroutine  --  here  labelled  MATMUL  --  will  compute  the 
value  of  the  product  of  two  matrices.  Separating  the  statement  that  one 
matrix  is  the  matrix  product  of  two  other  matrices  from  a  definition  of 
matrix  product  in  terms  of  a  subroutine  serves  to  separate  the  compu¬ 
tational  problem  of  determining  a  value  from  the  operational  problem 
of  deciding  what  to  do;  in  particular,  if  we  never  ask  for  the  value  of 
some  matrix  which  is  the  product  of  two  matrices,  we  never  need  to 
deduce  the  fact  that  the  subroutine  MATMUL  will  give  that  value. 

3  2 

The  transformation  from  a  polynomial  basis  [u  u  u  1]  to  the  end¬ 
point  derivative  basis  [F^F^G^G^]  is  given  simply  by  the  matrix  we 
have  called  M  and  is  represented  by  Figure  A. 2. 7.  The  transformation 
from  the  polynomial  basis  to  the  particular  basis  needed  to  generate  a 
surface  by  finite  difference  methods  is  a  function  of  the  number  and 
direction  of  the  curves  to  be  displayed  on  the  surface  and  the  number  of 
points  on  each  curve.  We  assume  a  subroutine  COMPUTEDISPLAYBASIS 
to  compute  the  value  of  the  matrix  for  this  transformation,  using  the 
material  of  subsections  2.2  and  3.5.  Figure  A. 2. 8  illustrates  these  con¬ 
cepts  and  would  be  applied  to  give  the  finite  difference  representation  of 
a  particular  surface  once  the  particular  family  of  parametric  curves  was 
chosen. 

In  a  similar  fashion  we  can  use  the  proposed  data  structure  to 
represent  the  algorithms  to  be  used  to  build  a  surface  specified  by  geo¬ 
metric  conditions.  Let  us  suppose  the  designer  wishes  to  define  a 
portion  of  a  cylinder  as  a  particular  kind  of  surface.  (We  speak  about  a 
portion  of  a  cylinder  since  we  cannot  draw  a  complete  cylinder  at  one 
surface  patch,  this  because  we  cannot  draw  a  complete  circle,  but  rather 
only  any  part  less  than  a  complete  circle.)  Granting  for  the  moment  that 
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such  a  definition  can  be  given.  Figure  A. 2.  9  represents  the  particular 
instance  GEORGE  of  a  cylinder,  including  the  fact  that  the  cylinder  is 
to  be  displayed  using  5  curves  in  the  u  =  constant  family  and  2  curves 
in  the  v  =  constant  family,  17  points  per  curve  in  the  former  and  33  points 
per  curve  in  the  latter.  The  particular  portion  of  a  cylinder  specified  is 
to  have  a  height  of  one  and  is  to  subtend  an  angle  of  270°.  In  order  to 
display  this  surface  we  need  to  go  from  the  specification  of  this  particu¬ 
lar  cylinder  to  a  tensor  representation  in  the  polynomial  basis;  by  the 
algorithm  expressed  in  Figure  A. 2. 8  we  can  then  find  the  particular  finite 
difference  representation  to  draw  the  surface. 

At  this  point  we  could  take  either  of  two  approaches.  The  first 
would  be  to  define  a  subroutine  MAKECYL  that  would  compute  the 
tensor  representation  of  a  cylinder  given  the  specified  dimensions  as 
arguments;  this  is  represented  in  Figure  A.  2. 10.  The  second  approach 
would  be  to  define  the  cylinder  as  a  particular  type  of  surface  of  revo¬ 
lution,  define  surface  of  revolution  as  a  particular  type  of  product 
surface,  and  then  define  a  product  surface  as  a  particular  type  of 
general  surface.  This  approach  has  the  advantage  of  directly  repre¬ 
senting  the  distinguishing  facts  about  the  cylinder  and  would  now  allow, 
for  instance,  a  sphere  to  be  readily  described  as  a  different  particular 
type  of  surface  of  revolution,  with  the  rest  of  the  machinery  directly 
applicable  to  its  generation. 

We  will  not  represent  the  characterization  of  the  cylinder  as  a 
surface  of  revolution  since  that  would  depend  on  some  unilluminating 
assumptions  about  the  representations  of  circle  arcs  and  straight  lines; 
we  will  indicate  how  the  rest  of  the  machinery  might  be  represented. 

Let  us  assume  a  surface  of  revolution  to  be  specified  by  a  generating 
curve  and  a  particular  arc  of  a  circle;  to  call  this  a  product  surface, 
we  need  only  relabel  the  generating  curve  and  circle  to  correspond  to 
the  two  curves  characterizing  a  product  surface  as  in  subsection  3.9  -- 
we  will  call  these  GENERATOR-ONE  and  GENERATOR-TWO.  This 
relabelling  is  indicated  by  the  representation  in  Figure  A. 2.11.  Then  to 
obtain  the  proper  general  tensor  representation  of  the  product  surface. 
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we  define  a  subroutine  PRODSURF  that  multiplies  the  matrices  repre¬ 
senting  the  curves  in  the  proper  way  to  give  the  value  of  the  tensor  -- 
notice  here  that  the  functional  basis  of  the  tensor  will  be  the  same  as 
the  single  basis  of  the  representations  of  the  two  generating  curves. 
These  concepts  are  represented  in  Figure  A  2  12,  completing  the  set  of 
deductions  necessary  to  draw  surfaces  of  revolution 

The  examples  just  discussed  are,  of  course,  not  complete  nor 
completely  accurate,  even  so  far  as  they  go.  In  particular,  the  change 
of  basis  transformations  would  be  more  accurately  expressed  in  some 
tensor-like  notation,  since  we  have  to  take  into  account  that  we  are 
dealing  with  two  bases,  the  functional  basis  of  the  curve  --  some  par¬ 
ticular  set  of  polynomials  --  and  the  basis  of  the  four-dimensional 
projective  space  in  which  the  object  is  represented.  This  latter  basis 
is  not  mentioned  at  all  in  the  examples,  but  must  of  course  be  included 
since  it  is  what  determines  the  position,  orientation,  and  perspective 
view  of  the  object.  For  instance,  we  would  like  to  allow  compound 
objects  composed  of  some  set  of  previously  defined  objects  in  some 
particular  orientation.  We  would  hope  that  a  spatial  transformation 
applied  to  such  a  compound  object  would  be  computed  by  applying  it 
first  to  the  transformations  expressing  the  internal  orientation  of  the 
individual  objects  in  the  compound  object,  and  then  applying  that  com¬ 
pound  transformation  to  the  representations  of  the  internal  objects. 

This  would,  of  course,  be  a  recursive  process,  since  a  component  of 
a  compound  object  might  be  a  compound  object  itself 

An  implementation  in  terms  of  this  data  structure  will  have  to 
provide  the  following  capabilities: 

1.  Representation  of  instances  of  relations  (ordered  triples). 

2.  Access  to  such  triples  using  any  one,  two,  or  all  three  of 
the  possible  elements  of  the  triple  as  addressing  criteria,  possibly 
producing  multiple  answers 

3.  Representations  of  numeric  or  other  values  of  objects. 

The  important  problem  to  be  solved  in  such  an  implementation  is 
that  all  relations  although  written  as  one-way  directed  arcs  are  to  be 
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two-way  arcs,  in  which  access  paths  exist  in  both  forward  and  backward 
directions.  Most  conventional  data  structures  do  not  provide  automati¬ 
cally  for  such  bi-directional  access  paths  in  a  symmetrical  manner. 

Even  fewer  permit  the  access  path  itself  to  be  a  possible  item  of  interest. 
The  LEAP  language  and  data  structure  [8,  18]  on  the  TX-2  at  the  M.  I.  T. 
Lincoln  Laboratory  is  an  example  of  a  fairly  efficient  way  to  represent 
this  type  of  information.  By  using  hashing  techniques,  it  stores  each 
instance  of  a  relation  in  three  ways  in  order  that  given  any  two  of  the 
components  of  a  triple,  there  is  a  virtually  direct  path  to  the  third  com¬ 
ponent. 

The  particular  class  of  problems  arising  in  the  context  of  surface 
design  involves  a  considerable  amount  of  purely  arithmetic  number 
manipulation,  such  as  in  multiplying  matrices  or  generating  displays. 

The  use  of  this  proposed  data  structure  and  the  overhead  of  extracting 
information  from  it  appear  likely  not  to  degrade  the  performance  of  the 
system  significantly,  and,  in  immensely  increasing  its  "cognitive" 
powers,  would  be  of  great  benefit. 
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Rational  parameteric  polynomial  functions  of  second  degree  or 
higher  provide  a  class  of  curves  including  all  conic  sections.  They 
can  be  generated  by  an  iterative  process  easily  implemented  in  soft¬ 
ware  or  hardware.  The  numerical  accuracy  of  the  process  is  analyzed. 
Algorithms  for  the  specification,  display,  and  modification  of  the  curve 
are  presented.  Such  curves  are  represented  in  a  homogeneous 
coordinate  formulation  convenient  for  computer  applications.  Three- 
dimensional  surfaces  composed  of  such  curves  are  similarly  conveni¬ 
ent  to  use.  Without  recourse  to  trigonometric  functions,  such  classi¬ 
cal  surfaces  as  spheres  and  toroids  can  be  readily  described.  The 
ease  with  which  translation,  rotation  and  projective  transformations 
can  be  applied  is  exhibited.  In  particular,  we  do  not  perform  such 
transformations  on  the  points  of  the  surface  to  be  displayed  --  upwards 
of  several  thousand  --  but  rather  upon  the  rather  small  set  of  numbers 
in  a  4  X  4  X  4  tensor  that  represents  the  surface.  These  surfaces  are 
intended  to  be  used  in  an  interactive,  freeform  computer-aided  design 
system.  In  this  direction  we  discuss  the  enforcing  of  continuity 
conditions  and  possible  data  structures  for  representing  the  surfaces. 
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